cloudy
trunk
|
00001 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and 00002 * others. For conditions of distribution and use see copyright notice in license.txt */ 00003 /*plot master routine to generate some sort of plot */ 00004 #include "cddefines.h" 00005 #include "iterations.h" 00006 #include "called.h" 00007 #define IHI 59 00008 #define IWID 121 00009 #include "input.h" 00010 #include "rfield.h" 00011 #include "trace.h" 00012 #include "radius.h" 00013 #include "geometry.h" 00014 #include "opacity.h" 00015 #include "dense.h" 00016 #include "hcmap.h" 00017 #include "plot.h" 00018 00019 /*pltcon generate plot of continuum array */ 00020 STATIC void pltcon(long int np, 00021 const char *chCall); 00022 00023 /*pltmap generate plot of heating and cooling map */ 00024 STATIC void pltmap(long int np, 00025 const char *chCall); 00026 00027 /*pltopc generate plot of local gas opacity */ 00028 STATIC void pltopc(long int np, 00029 const char *chCall); 00030 00031 /* this is the base routine that actually makes the plots, called by above */ 00032 STATIC void pltr(realnum[],realnum[],long,double,double,double,double, 00033 char,char*,long,bool); 00034 00035 00036 void plot(const char *chCall) 00037 { 00038 long int np; 00039 00040 DEBUG_ENTRY( "plot()" ); 00041 00042 /* return if this is not the last iteration, or a plot not required, 00043 * or we are not speaking */ 00044 if( !plotCom.lgPlotON || !called.lgTalk ) 00045 { 00046 return; 00047 } 00048 00049 if( !iterations.lgLastIt && (strcmp(chCall,"FIRST") != 0) ) 00050 { 00051 return; 00052 } 00053 00054 /* loop over all the requested plots */ 00055 for( np=0; np < plotCom.nplot; np++ ) 00056 { 00057 /* series of tests to determine which type of plot we will do */ 00058 if( strcmp(plotCom.chPType[np]," MAP") == 0 ) 00059 { 00060 /* thermal map */ 00061 pltmap(np,chCall); 00062 } 00063 else if( strcmp(plotCom.chPType[np] ,"CONT") == 0 || 00064 strcmp(plotCom.chPType[np] ,"CRAW") == 0 || 00065 strcmp(plotCom.chPType[np] ,"DIFF") == 0 || 00066 strcmp(plotCom.chPType[np] ,"REFL") == 0 || 00067 strcmp(plotCom.chPType[np] ,"EMIT") == 0 || 00068 strcmp(plotCom.chPType[np] ,"CPHT") == 0 || 00069 strcmp(plotCom.chPType[np] ,"OUTW") == 0 ) 00070 { 00071 /* this is a contiuum plot of some kind */ 00072 pltcon(np,chCall); 00073 } 00074 00075 else if( 00076 strcmp(plotCom.chPType[np] ,"OPAA") == 0 || 00077 strcmp(plotCom.chPType[np] ,"OPAS") == 0 || 00078 strcmp(plotCom.chPType[np] ,"OPAT") == 0 ) 00079 { 00080 /* absorption, scattering, or total opacity */ 00081 pltopc(np,chCall); 00082 } 00083 else 00084 { 00085 fprintf( ioQQQ, " PLOT type=%4.4s not known. STOP\n", 00086 plotCom.chPType[np] ); 00087 cdEXIT(EXIT_FAILURE); 00088 } 00089 } 00090 00091 return; 00092 } 00093 00094 /*pltcon generate plot of continuum array */ 00095 00096 STATIC void pltcon( 00097 long int np, 00098 const char *chCall) 00099 { 00100 char chSymPlt2[3], 00101 chXtitle[23]; 00102 char chSym, 00103 chSymPlt1; 00104 long int i; 00105 double contin, 00106 ymin2; 00107 static double xmax, 00108 xmin, 00109 ymax, 00110 ymin; 00111 static realnum *y/*[rfield.nupper]*/, 00112 *y2/*[rfield.nupper]*/; 00113 00114 DEBUG_ENTRY( "pltcon()" ); 00115 00116 if( strcmp(chCall,"FIRST") == 0 ) 00117 { 00118 return; 00119 } 00120 00121 xmin = rfield.anulog[0]; 00122 xmin = MAX2((double)plotCom.pltxmn[np],xmin); 00123 xmax = rfield.anulog[rfield.nflux-1]; 00124 xmax = MIN2(xmax,(double)plotCom.pltxmx[np]); 00125 00126 if( plotCom.lgPltTrace[np] ) 00127 { 00128 fprintf( ioQQQ, " XMIN, XMAX=%12.4e%12.4e NFLUX=%4ld\n", 00129 xmin, xmax, rfield.nflux ); 00130 } 00131 00132 if( strcmp(plotCom.chPType[np],"REFL") == 0 && geometry.lgSphere ) 00133 { 00134 fprintf( ioQQQ, " Reflected continuum not computed when SPHERE set.\n" ); 00135 return; 00136 } 00137 00138 y = (realnum*)MALLOC((size_t)rfield.nupper*sizeof(realnum) ); 00139 y2 = (realnum*)MALLOC((size_t)rfield.nupper*sizeof(realnum) ); 00140 00141 /* these will be the default symbols for first and second plot */ 00142 chSymPlt1 = '.'; 00143 strcpy( chSymPlt2, "o " ); 00144 ymin = FLT_MAX; 00145 ymin2 = FLT_MAX; 00146 ymax = -FLT_MAX; 00147 for( i=0; i < rfield.nflux; i++ ) 00148 { 00149 if( (double)rfield.anulog[i] > xmin && (double)rfield.anulog[i] < xmax ) 00150 { 00151 if( strcmp(plotCom.chPType[np],"CONT") == 0 ) 00152 { 00153 y[i] = (realnum)log10(MAX2(rfield.flux_total_incident[i]/rfield.widflx[i]* 00154 rfield.anu2[i],1e-37)); 00155 /* >>chng 01 jul 13, add rfield.ConEmitReflec[i] */ 00156 contin = rfield.flux[i] + rfield.ConEmitOut[i]*geometry.covgeo + rfield.ConEmitReflec[i]; 00157 y2[i] = (realnum)MAX2((contin/rfield.widflx[i]+(rfield.outlin[i]+rfield.outlin_noplot[i])/ 00158 rfield.anu[i]*geometry.covgeo)* 00159 rfield.anu2[i]*radius.r1r0sq,1e-37); 00160 y2[i] = (realnum)log10(y2[i]); 00161 } 00162 else if( strcmp(plotCom.chPType[np],"CPHT") == 0 ) 00163 { 00164 /* plot continuum as photons */ 00165 y[i] = (realnum)log10(MAX2(rfield.flux_total_incident[i]/rfield.widflx[i], 00166 1e-37)); 00167 contin = rfield.flux[i] + rfield.ConEmitOut[i]*geometry.covgeo/rfield.widflx[i]; 00168 y2[i] = (realnum)MAX2((contin+(rfield.outlin[i]+rfield.outlin[i])/rfield.anu[i]* 00169 geometry.covgeo)* radius.r1r0sq,1e-37); 00170 y2[i] = (realnum)log10(y2[i]); 00171 } 00172 else if( strcmp(plotCom.chPType[np],"REFL") == 0 ) 00173 { 00174 /* plot "reflected" continuum from last zone only */ 00175 y[i] = (realnum)log10(MAX2((rfield.ConEmitReflec[i]/rfield.widflx[i]+ 00176 rfield.reflin[i])*rfield.anu2[i],1e-37)); 00177 y2[i] = y[i]; 00178 } 00179 else if( strcmp(plotCom.chPType[np],"EMIT") == 0 ) 00180 { 00181 /* plot "emitted" continuum from both sides of cloud */ 00182 y[i] = (realnum)log10(MAX2( 00183 ((rfield.ConEmitReflec[i]+rfield.ConEmitOut[i])/ 00184 rfield.widflx[i]+ 00185 (rfield.outlin[i]+rfield.outlin_noplot[i]+rfield.reflin[i])/rfield.anu[i] )* 00186 rfield.anu2[i],1e-37)); 00187 y2[i] = y[i]; 00188 } 00189 else if( strcmp(plotCom.chPType[np],"OUTW") == 0 ) 00190 { 00191 /* plot outward and attenuated incident continuum */ 00192 chSymPlt1 = 'i'; 00193 y[i] = (realnum)log10(MAX2(rfield.flux[i]*opac.opacity_abs[i], 00194 1e-37)); 00195 strcpy( chSymPlt2, "o " ); 00196 y2[i] = (realnum)log10(MAX2((rfield.outlin[i]+rfield.outlin_noplot[i]+rfield.ConEmitOut[i])* 00197 opac.opacity_abs[i],1e-37)); 00198 } 00199 else if( strcmp(plotCom.chPType[np],"DIFF") == 0 ) 00200 { 00201 /* plot "diffuse" continuum from last zone only */ 00202 y[i] = (realnum)log10(MAX2(rfield.ConEmitLocal[nzone][i]*rfield.anu2[i]/ 00203 rfield.widflx[i],1e-37)); 00204 y2[i] = y[i]; 00205 } 00206 else if( strcmp(plotCom.chPType[np],"CRAW") == 0 ) 00207 { 00208 y[i] = (realnum)log10(MAX2(rfield.flux_total_incident[i],1e-37)); 00209 y2[i] = (realnum)MAX2((rfield.flux[i]+ 00210 rfield.otscon[i]+rfield.otslin[i]+rfield.outlin[i]+rfield.outlin_noplot[i]+ 00211 rfield.ConEmitOut[i])*radius.r1r0sq,1e-37); 00212 y2[i] = (realnum)log10(y2[i]); 00213 } 00214 00215 if( y[i] > -36.9 ) 00216 { 00217 ymin = MIN2(ymin,(double)y[i]); 00218 } 00219 00220 if( y2[i] > -36.9 ) 00221 { 00222 ymin2 = MIN2(ymin2,(double)y2[i]); 00223 } 00224 00225 ymax = MAX2(ymax,(double)y[i]); 00226 ymax = MAX2(ymax,(double)y2[i]); 00227 } 00228 } 00229 00230 if( trace.lgTrace ) 00231 { 00232 fprintf( ioQQQ, " PLOT called for the first time, YMAX, MIN=%10.2e%10.2e\n", 00233 ymax, ymin ); 00234 } 00235 00236 /* lower min by at most 5 dex below peak */ 00237 ymin2 = MAX3(ymax-5.,-35.,ymin2); 00238 00239 /* make sure there is room at the bottom */ 00240 ymin = MIN3(ymin2,ymin,ymax-1.); 00241 00242 /* emitted continuum is thermal, so goes to zero */ 00243 if( strcmp(plotCom.chPType[np],"EMIT") == 0 ) 00244 { 00245 ymin = MAX2(ymin,ymax-4.); 00246 } 00247 00248 if( plotCom.lgPltTrace[np] ) 00249 { 00250 fprintf( ioQQQ, " YMAX, MIN=%14.4e%14.4e Npnts=%4ld\n", 00251 ymax 00252 , ymin, rfield.nflux ); 00253 } 00254 strcpy( chXtitle, "Log(nu fnu) vs LOG(nu)" ); 00255 00256 chSym = chSymPlt1; 00257 00258 pltr(rfield.anulog,y,rfield.nflux,xmin,xmax,ymin,ymax,chSym,chXtitle 00259 ,1,plotCom.lgPltTrace[np]); 00260 00261 chSym = chSymPlt2[0]; 00262 00263 pltr(rfield.anulog,y2,rfield.nflux,xmin,xmax,ymin,ymax,chSym,chXtitle 00264 ,3,plotCom.lgPltTrace[np]); 00265 00266 free( y ); 00267 free( y2 ); 00268 return; 00269 } 00270 00271 /*pltmap generate plot of heating and cooling map */ 00272 00273 STATIC void pltmap( 00274 long int np, 00275 const char *chCall) 00276 { 00277 char chXtitle[23]; 00278 static bool lgTlkSav; 00279 char chSym; 00280 00281 long int i; 00282 00283 static double xmax, 00284 xmin, 00285 ymax, 00286 ymin; 00287 00288 DEBUG_ENTRY( "pltmap()" ); 00289 00290 if( strcmp(chCall,"FIRST") == 0 ) 00291 { 00292 return; 00293 } 00294 00295 lgTlkSav = called.lgTalk; 00296 called.lgTalk = false; 00297 hcmap.lgMapBeingDone = true; 00298 hcmap.RangeMap[0] = (realnum)pow((realnum)10.f,plotCom.pltxmn[np]); 00299 hcmap.RangeMap[1] = (realnum)pow((realnum)10.f,plotCom.pltxmx[np]); 00300 map_do(ioQQQ, " map"); 00301 called.lgTalk = lgTlkSav; 00302 00303 for( i=0; i < hcmap.nmap; i++ ) 00304 { 00305 hcmap.temap[i] = (realnum)log10(hcmap.temap[i]); 00306 } 00307 00308 xmin = MIN2(hcmap.temap[0],hcmap.temap[hcmap.nmap-1]); 00309 xmin = MAX2((double)plotCom.pltxmn[np],xmin); 00310 xmax = MAX2(hcmap.temap[0],hcmap.temap[hcmap.nmap-1]); 00311 xmax = MIN2(xmax,(double)plotCom.pltxmx[np]); 00312 00313 if( plotCom.lgPltTrace[np] ) 00314 { 00315 fprintf( ioQQQ, " xmin, xmax=%12.4e%12.4e nmap=%4ld\n", 00316 xmin, xmax, hcmap.nmap ); 00317 } 00318 00319 ymin = FLT_MAX; 00320 ymax = -FLT_MAX; 00321 00322 for( i=0; i < hcmap.nmap; i++ ) 00323 { 00324 if( (double)hcmap.temap[i] > xmin && (double)hcmap.temap[i] < xmax ) 00325 { 00326 hcmap.hmap[i] = (realnum)log10(MAX2(hcmap.hmap[i],1e-35)); 00327 hcmap.cmap[i] = (realnum)log10(MAX2(hcmap.cmap[i],1e-35)); 00328 if( hcmap.cmap[i] > -34. ) 00329 { 00330 ymin = MIN3(ymin,hcmap.hmap[i],hcmap.cmap[i]); 00331 } 00332 else 00333 { 00334 ymin = MIN2(ymin,(double)hcmap.hmap[i]); 00335 } 00336 ymax = MAX3(ymax,hcmap.hmap[i],hcmap.cmap[i]); 00337 } 00338 } 00339 00340 if( trace.lgTrace ) 00341 { 00342 fprintf( ioQQQ, " PLOT called for the first time, YMAX, MIN=%10.2e%10.2e\n", 00343 ymax, ymin ); 00344 } 00345 00346 if( plotCom.lgPltTrace[np] ) 00347 { 00348 fprintf( ioQQQ, " YMAX, MIN=%14.4e%14.4e Npnts=%4ld\n", 00349 ymax, ymin, hcmap.nmap ); 00350 } 00351 00352 chSym = 'H'; 00353 strcpy( chXtitle, "heating - cooling v te" ); 00354 00355 pltr(hcmap.temap,hcmap.hmap,hcmap.nmap,xmin,xmax,ymin,ymax,chSym, 00356 chXtitle,1,plotCom.lgPltTrace[np]); 00357 00358 chSym = 'C'; 00359 00360 pltr(hcmap.temap,hcmap.cmap,hcmap.nmap,xmin,xmax,ymin,ymax,chSym, 00361 chXtitle,3,plotCom.lgPltTrace[np]); 00362 00363 return; 00364 } 00365 00366 /*pltopc generate plot of local gas opacity */ 00367 STATIC void pltopc( 00368 long int np, 00369 const char *chCall) 00370 { 00371 char chXtitle[23]; 00372 char chSym; 00373 long int i; 00374 double arg1, 00375 arg2; 00376 static double xmax, 00377 xmin, 00378 ymax, 00379 ymin; 00380 static realnum *y/*[rfield.nupper]*/, 00381 *y2/*[rfield.nupper]*/; 00382 00383 DEBUG_ENTRY( "pltopc()" ); 00384 00385 if( strcmp(chCall,"FIRST") == 0 ) 00386 { 00387 return; 00388 } 00389 00390 y = (realnum*)MALLOC((size_t)rfield.nupper*sizeof(realnum) ); 00391 y2 = (realnum*)MALLOC((size_t)rfield.nupper*sizeof(realnum) ); 00392 00393 xmin = rfield.anulog[0]; 00394 xmin = MAX2((double)plotCom.pltxmn[np],xmin); 00395 xmax = rfield.anulog[rfield.nflux-1]; 00396 xmax = MIN2(xmax,(double)plotCom.pltxmx[np]); 00397 00398 if( plotCom.lgPltTrace[np] ) 00399 { 00400 fprintf( ioQQQ, " XMIN, XMAX=%12.4e%12.4e NFLUX=%4ld\n", 00401 xmin, xmax, rfield.nflux ); 00402 } 00403 00404 ymin = FLT_MAX; 00405 ymax = -FLT_MAX; 00406 00407 for( i=0; i < rfield.nflux; i++ ) 00408 { 00409 if( strcmp(plotCom.chPType[np],"OPAA") == 0 ) 00410 { 00411 /* absorption opacity */ 00412 arg1 = opac.opacity_abs_savzon1[i]; 00413 arg2 = opac.opacity_abs[i]; 00414 } 00415 00416 else if( strcmp(plotCom.chPType[np],"OPAS") == 0 ) 00417 { 00418 /* scattering opacity */ 00419 arg1 = opac.opacity_sct_savzon1[i]; 00420 arg2 = opac.opacity_sct[i]; 00421 } 00422 00423 else if( strcmp(plotCom.chPType[np],"OPAT") == 0 ) 00424 { 00425 /* total opacity */ 00426 arg1 = opac.opacity_abs_savzon1[i] + opac.opacity_sct_savzon1[i]; 00427 arg2 = opac.opacity_abs[i] + opac.opacity_sct[i]; 00428 } 00429 00430 else 00431 { 00432 /* this cannot happen since type was set to one of above */ 00433 fprintf( ioQQQ, " pltopc type=%4.4s not known. STOP\n", 00434 plotCom.chPType[np] ); 00435 cdEXIT(EXIT_FAILURE); 00436 } 00437 00438 y[i] = (realnum)log10(MAX2(arg1/dense.gas_phase[ipHYDROGEN],1e-35)); 00439 y2[i] = (realnum)log10(MAX2(arg2/dense.gas_phase[ipHYDROGEN],1e-35)); 00440 00441 if( (double)rfield.anulog[i] > xmin && (double)rfield.anulog[i] < xmax ) 00442 { 00443 ymin = MIN3(ymin,y[i],y2[i]); 00444 ymax = MAX3(ymax,y[i],y2[i]); 00445 } 00446 } 00447 00448 if( trace.lgTrace ) 00449 { 00450 fprintf( ioQQQ, " PLOT called for the first time, YMAX, MIN=%10.2e%10.2e\n", 00451 ymax, ymin ); 00452 } 00453 00454 /* lower min by factor of 10 to show absorption in next plot */ 00455 ymin = MAX2(ymin-1.,-35.); 00456 ymax += 1.; 00457 if( plotCom.lgPltTrace[np] ) 00458 { 00459 fprintf( ioQQQ, " YMAX, MIN=%14.4e%14.4e Npnts=%4ld\n", 00460 ymax, ymin, rfield.nflux ); 00461 } 00462 00463 strcpy( chXtitle, "Log(opacity) vs log(n)" ); 00464 00465 chSym = '.'; 00466 pltr(rfield.anulog,y,rfield.nflux,xmin,xmax,ymin,ymax,chSym,chXtitle 00467 ,1,plotCom.lgPltTrace[np]); 00468 00469 chSym = 'o'; 00470 pltr(rfield.anulog,y2,rfield.nflux,xmin,xmax,ymin,ymax,chSym,chXtitle 00471 ,3,plotCom.lgPltTrace[np]); 00472 00473 free(y); 00474 free(y2); 00475 return; 00476 } 00477 00478 /*pltr core plotting routine for generating line printer plots */ 00479 00480 STATIC void pltr( 00481 /* the x-axis */ 00482 realnum x[], 00483 /* the y-axi */ 00484 realnum y[], 00485 /* number of points */ 00486 long int npnts, 00487 /* mins and maxs, log of min and max of x-axis */ 00488 double xmin, 00489 double xmax, 00490 double ymin, 00491 double ymax, 00492 /* plot symbol */ 00493 char chSym, 00494 char *chXtitle, 00495 long int itim, 00496 bool lgTrace) 00497 { 00498 static char chPage[59][122]; 00499 00500 long int i, 00501 ix, 00502 iy, 00503 j, 00504 nc; 00505 00506 /* the max number of decades we can plot */ 00507 # define NDECAD 18 00508 00509 static long int jpnt[NDECAD], 00510 lowx, 00511 lx; 00512 00513 static double xdec, 00514 xinc, 00515 ydown, 00516 yinc; 00517 00518 /* this is log of smallestnumer in following set */ 00519 const realnum xAxisMin = -8.f; 00520 00521 static char chLab[NDECAD][5]={"1E-8","1E-7","1E-6","1E-5", 00522 "1E-4",".001","0.01"," 0.1"," 1 ", 00523 " 10 "," 100","1000","1E4 ","1E5 ","1E6 ","1E7 ","1E8 ","1E9 "}; 00524 00525 DEBUG_ENTRY( "pltr()" ); 00526 00527 /* ITIM=1, first call, =2 intermediate calls, =3 for last call*/ 00528 if( itim == 1 ) 00529 { 00530 /* first call, set left border of plot and clear out array */ 00531 for( i=1; i < IHI; i++ ) 00532 { 00533 chPage[i][0] = 'l'; 00534 for( j=1; j < IWID; j++ ) 00535 { 00536 chPage[i][j] = ' '; 00537 } 00538 } 00539 00540 /* centered label for plot */ 00541 strcpy( chPage[1], " " ); 00542 strcat( chPage[1], chXtitle ); 00543 strcat( chPage[1], input.chTitle ); 00544 00545 /* one dex increments in x and y marked special */ 00546 i = 1; 00547 ydown = 0.; 00548 yinc = (realnum)(IHI-2)/(ymax - ymin); 00549 nc = 0; 00550 00551 while( i <= IHI && nc < 200 ) 00552 { 00553 chPage[i-1][1] = '-'; 00554 ydown += yinc; 00555 i = (long)(ydown + 1); 00556 nc += 1; 00557 } 00558 00559 /* bottom increments of plot */ 00560 for( i=0; i < IWID; i++ ) 00561 { 00562 chPage[IHI-1][i] = '-'; 00563 } 00564 00565 if( xmin < xAxisMin ) 00566 { 00567 fprintf(ioQQQ," plts: xmin is less than min value in array\n"); 00568 cdEXIT(EXIT_FAILURE); 00569 } 00570 /* LX is pointer to label for number in x-axis in chLab */ 00571 if( xmin < 0. ) 00572 { 00573 lx = (long)(4.999-fabs(xmin)); 00574 /* lx is the offset within the array of x-axis values */ 00575 /* >>chng 99 jun 11 change to allow any min value of x-axis */ 00576 lx = (long)(fabs(xAxisMin)-0.001-fabs(xmin)); 00577 lx = MAX2(0,lx); 00578 /* this is lowest decade on plot */ 00579 xdec = -floor(fabs(xmin)+1e-5); 00580 } 00581 else 00582 { 00583 double aa; 00584 lx = (long)MAX2(0.,4.+xmin); 00585 /* lx is the offset within the array of x-axis values */ 00586 lx = (long)MAX2(0.,4.+xmin); 00587 /* >>chng 99 jun 11 change to allow any min value of x-axis */ 00588 aa = fabs(xAxisMin); 00589 lx = (long)MAX2(0., aa-1. + xmin ); 00590 xdec = floor(xmin+1e-5); 00591 } 00592 00593 lowx = lx + 1; 00594 xinc = (realnum)(IWID-1)/(xmax - xmin); 00595 i = (long)MAX2(1.,(xdec-xmin)*xinc+1.); 00596 nc = 0; 00597 00598 while( i < IWID && nc < 100 ) 00599 { 00600 chPage[IHI-2][i - 1] = 'l'; 00601 00602 /* fix position of x labels */ 00603 lx = MIN2(lx+1,NDECAD); 00604 00605 /* slight offset to center label */ 00606 jpnt[lx-1] = MAX2(0,i-3); 00607 jpnt[lx-1] = MIN2((long)IWID-4,jpnt[lx-1]); 00608 xdec += 1.; 00609 i = (long)MAX2(1.,(xdec-xmin)*xinc+1.); 00610 nc += 1; 00611 } 00612 } 00613 00614 /* everything falls down through here */ 00615 /* now fill in data, symbol is chSym */ 00616 for( i=0; i < npnts; i++ ) 00617 { 00618 if( (double)x[i] > xmin && (double)x[i] < xmax ) 00619 { 00620 iy = (long)(IHI - MAX2(y[i]-ymin,0.)*yinc); 00621 iy = MAX2(1,iy); 00622 ix = (long)((x[i] - xmin)*xinc + 1); 00623 00624 if( lgTrace ) 00625 { 00626 fprintf( ioQQQ, " x, y, ix, iy=%7.3f%7.3f%4ld%4ld\n", 00627 x[i], y[i], ix, iy ); 00628 } 00629 chPage[iy-1][ix - 1] = chSym; 00630 } 00631 } 00632 00633 if( itim == 3 ) 00634 { 00635 /* make the plot */ 00636 fprintf( ioQQQ, "1\n" ); 00637 for( i=1; i < IHI; i++ ) 00638 { 00639 fprintf( ioQQQ, " %121.121s\n", chPage[i] ); 00640 } 00641 00642 /* now put on label for X-axis */ 00643 for( i=0; i < IWID; i++ ) 00644 { 00645 chPage[0][i] = ' '; 00646 } 00647 00648 for( i=lowx-1; i < lx; i++ ) 00649 { 00650 /* copy the four char of the numeric string */ 00651 strncpy(chPage[0]+jpnt[i] , chLab[i+1] , 4); 00652 } 00653 fprintf( ioQQQ, " %121.121s\n", chPage[0] ); 00654 } 00655 return; 00656 }