00001
00002
00003
00004
00005
00006 #include "cddefines.h"
00007 #include "input.h"
00008 #include "conv.h"
00009 #include "optimize.h"
00010 #include "iso.h"
00011 #include "called.h"
00012 #include "atmdat.h"
00013 #include "hcmap.h"
00014 #include "thermal.h"
00015 #include "pressure.h"
00016 #include "struc.h"
00017 #include "wind.h"
00018 #include "h2.h"
00019 #include "colden.h"
00020 #include "dense.h"
00021 #include "lines.h"
00022 #include "secondaries.h"
00023 #include "radius.h"
00024 #include "version.h"
00025 #include "hmi.h"
00026 #include "prt.h"
00027 #include "grainvar.h"
00028 #include "atomfeii.h"
00029 #include "cddrive.h"
00030 #include "elementnames.h"
00031 #include "monitor_results.h"
00032 #include "taulines.h"
00033 #include "grid.h"
00034 #include "lines_service.h"
00035 #include "parser.h"
00036
00037 bool lgMonitorsOK , lgBigBotch, lgPrtSciNot;
00038 t_monitorresults MonitorResults;
00039
00040
00041 static bool lgInitDone=false ,
00042
00043 lgSpaceAllocated=false;
00044
00045
00046 static const int NASSERTS = 200;
00047
00048
00049 static realnum ErrorDefault;
00050
00051
00052 static realnum ErrorDefaultPerformance;
00053
00054
00055 static const int NCHAR = 5;
00056 static char **chAssertLineLabel;
00057
00058
00059 static char *chAssertLimit;
00060
00061
00062 static char **chAssertType;
00063
00064
00065 static double *AssertQuantity ,*AssertQuantity2 ,*AssertError, **Param;
00066
00067
00068 static int *lgQuantityLog;
00069 static long nAsserts=0;
00070 static realnum *wavelength;
00071
00072 inline double ForcePass( char chAssertLimit1 )
00073 {
00074
00075 if( chAssertLimit1 == '=' )
00076 return 0.;
00077 else if( chAssertLimit1 == '<' )
00078 return 1.;
00079 else if( chAssertLimit1 == '>' )
00080 return -1.;
00081 else
00082 TotalInsanity();
00083 }
00084
00085
00086
00087 void InitMonitorResults(void)
00088 {
00089
00090
00091
00092 lgInitDone = true;
00093 nAsserts = 0;
00094
00095
00096
00097 ErrorDefault = 0.05f;
00098
00099 ErrorDefaultPerformance = 0.2f;
00100 }
00101
00102
00103
00104 void ParseMonitorResults(Parser &p)
00105 {
00106 long i ,
00107 nelem,
00108 n2;
00109 char chLabel[INPUT_LINE_LENGTH];
00110
00111 DEBUG_ENTRY( "ParseMonitorResults()" );
00112
00113 if( !lgInitDone )
00114 {
00115 fprintf( ioQQQ, " ParseMonitorResults called before InitAsserResults\n" );
00116 cdEXIT(EXIT_FAILURE);
00117 }
00118
00119
00120 if( !lgSpaceAllocated )
00121 {
00122
00123
00124 lgSpaceAllocated = true;
00125
00126
00127 chAssertLineLabel = ((char **)MALLOC(NASSERTS*sizeof(char *)));
00128
00129
00130 chAssertType = ((char **)MALLOC(NASSERTS*sizeof(char *)));
00131
00132
00133 Param = ((double **)MALLOC(NASSERTS*sizeof(double * )));
00134
00135
00136 for( i=0; i<NASSERTS; ++i )
00137 {
00138 chAssertLineLabel[i] = ((char *)MALLOC(NCHAR*sizeof(char )));
00139 strcpy(chAssertLineLabel[i],"unkn" );
00140
00141
00142 chAssertType[i] = ((char *)MALLOC(3*sizeof(char )));
00143 strcpy(chAssertType[i],"un" );
00144
00145
00146 Param[i] = ((double *)MALLOC(5*sizeof(double )));
00147 }
00148
00149
00150 AssertQuantity = ((double *)MALLOC(NASSERTS*sizeof(double )));
00151
00152 AssertQuantity2 = ((double *)MALLOC(NASSERTS*sizeof(double )));
00153
00154
00155 AssertError = ((double *)MALLOC(NASSERTS*sizeof(double )));
00156
00157
00158 wavelength = ((realnum *)MALLOC(NASSERTS*sizeof(realnum )));
00159
00160
00161 lgQuantityLog = ((int *)MALLOC(NASSERTS*sizeof(int )));
00162
00163
00164 chAssertLimit = ((char *)MALLOC(NASSERTS*sizeof(char )));
00165 }
00166
00167
00168
00169 vector<double> AssertVector;
00170 if( p.nMatch( "GRID" ) )
00171 {
00172 ASSERT( optimize.nOptimiz >= 0 );
00173 AssertVector.resize( optimize.nOptimiz+1 );
00174
00175
00176 p.GetQuote( chLabel, true );
00177
00178 bool lgEOF;
00179 input_readvector( chLabel, get_ptr(AssertVector), optimize.nOptimiz+1, &lgEOF );
00180 if( lgEOF )
00181 fprintf(ioQQQ,"PROBLEM the file %s does not have enough values. sorry\n", chLabel );
00182 }
00183
00184
00185 strcpy( chAssertType[nAsserts] , "un" );
00186
00187
00188
00189 lgQuantityLog[nAsserts] = false;
00190
00191
00192 if( p.nMatch("<" ) )
00193 {
00194 chAssertLimit[nAsserts] = '<';
00195 }
00196 else if( p.nMatch(">" ) )
00197 {
00198 chAssertLimit[nAsserts] = '>';
00199 }
00200 else
00201 {
00202 chAssertLimit[nAsserts] = '=';
00203 }
00204
00205
00206
00207 if( p.nMatch(" SET" ) )
00208 {
00209
00210
00211
00212
00213
00214 if( nAsserts >0 )
00215 fprintf(ioQQQ," The default monitor error is being changed after"
00216 " some asserts were entered. \n This will only affect asserts"
00217 " that come after this command.\n");
00218 --nAsserts;
00219
00220 if( p.nMatch("ERRO" ) )
00221 {
00222 if( p.nMatch("PERF" ) )
00223 {
00224
00225 ErrorDefaultPerformance =
00226 (realnum)p.FFmtRead();
00227 if( p.lgEOL() )
00228 p.NoNumb("error");
00229 }
00230
00231 ErrorDefault =
00232 (realnum)p.FFmtRead();
00233 if( p.lgEOL() )
00234 p.NoNumb("error");
00235 }
00236 else
00237 {
00238
00239 fprintf( ioQQQ,
00240 " I could not identify an option on this ASSERT SET XXX command.\n");
00241 fprintf( ioQQQ, " Sorry.\n" );
00242 cdEXIT(EXIT_FAILURE);
00243 }
00244 }
00245
00246
00247 else if( p.nMatch("IONI" ) )
00248 {
00249
00250
00251
00252
00253
00254
00255 if( p.nMatch("VOLU" ) )
00256 {
00257 strcpy( chAssertType[nAsserts] , "F " );
00258 }
00259 else
00260 {
00261
00262 strcpy( chAssertType[nAsserts] , "f " );
00263 }
00264
00265
00266 if( (nelem = p.GetElem()) < 0 )
00267 {
00268 fprintf( ioQQQ,
00269 " I could not identify an element name on this line.\n");
00270 fprintf( ioQQQ, " Sorry.\n" );
00271 cdEXIT(EXIT_FAILURE);
00272 }
00273 ASSERT( nelem>= 0);
00274 ASSERT( nAsserts>= 0);
00275
00276
00277 strcpy( chAssertLineLabel[nAsserts], elementnames.chElementNameShort[nelem] );
00278
00279
00280 wavelength[nAsserts] = (realnum)p.FFmtRead();
00281 if( p.lgEOL() )
00282 {
00283 p.NoNumb("ionization stage");
00284 }
00285
00286 if( wavelength[nAsserts] < 1 || wavelength[nAsserts] > nelem+2 )
00287 {
00288 fprintf( ioQQQ,
00289 " The ionization stage is inappropriate for this element.\n");
00290 fprintf( ioQQQ, " Sorry.\n" );
00291 cdEXIT(EXIT_FAILURE);
00292 }
00293
00294 if( p.nMatch( "GRID" ) )
00295 {
00296 AssertQuantity[nAsserts] = AssertVector[optimize.nOptimiz];
00297 }
00298 else
00299 {
00300
00301
00302 AssertQuantity[nAsserts] = p.FFmtRead();
00303 if( p.lgEOL() )
00304 p.NoNumb("ionization fraction");
00305 }
00306
00307
00308 AssertError[nAsserts] = p.FFmtRead();
00309 if( p.lgEOL() )
00310 AssertError[nAsserts] = ErrorDefault;
00311
00312
00313
00314 if( AssertQuantity[nAsserts] <= 0. )
00315 {
00316
00317 AssertQuantity[nAsserts] =
00318 pow(10.,AssertQuantity[nAsserts] );
00319
00320 lgQuantityLog[nAsserts] = true;
00321 }
00322 else if( AssertQuantity[nAsserts] > 1. )
00323 {
00324 fprintf( ioQQQ,
00325 " The ionization fraction must be less than one.\n");
00326 fprintf( ioQQQ, " Sorry.\n" );
00327 cdEXIT(EXIT_FAILURE);
00328 }
00329
00330
00331 if( fabs(AssertQuantity[nAsserts]) <= SMALLDOUBLE )
00332 {
00333 fprintf( ioQQQ,
00334 " The ionization ionization fraction is too small, or zero. Check input\n");
00335 fprintf( ioQQQ, " Sorry.\n" );
00336 cdEXIT(EXIT_FAILURE);
00337 }
00338 }
00339
00340
00341 else if( p.nMatch("MOLE" )&&p.nMatch("FRAC" ) )
00342 {
00343
00344
00345
00346
00347
00348
00350 if( p.nMatch("VOLU" ) )
00351 {
00352 strcpy( chAssertType[nAsserts] , "MF" );
00353 }
00354 else
00355 {
00356
00357 strcpy( chAssertType[nAsserts] , "mf" );
00358 }
00359
00360 if( p.nMatchErase(" H2 " ) )
00361 {
00362 strcpy( chAssertLineLabel[nAsserts], "H2 " );
00363
00364 }
00365 else if( p.nMatch(" CO " ) )
00366 {
00367 strcpy( chAssertLineLabel[nAsserts], "CO " );
00368 }
00369 else
00370 {
00371 fprintf( ioQQQ,
00372 " I could not identify CO or H2 on this line.\n");
00373 fprintf( ioQQQ, " Sorry.\n" );
00374 cdEXIT(EXIT_FAILURE);
00375 }
00376
00377
00378 wavelength[nAsserts] = 0;
00379
00380
00381 AssertQuantity[nAsserts] = p.FFmtRead();
00382 if( p.lgEOL() )
00383 {
00384 p.NoNumb("molecular fraction");
00385 }
00386 if( AssertQuantity[nAsserts] <= 0. )
00387 {
00388
00389 AssertQuantity[nAsserts] =
00390 pow(10.,AssertQuantity[nAsserts] );
00391 }
00392
00393
00394
00395 AssertError[nAsserts] = p.FFmtRead();
00396 if( p.lgEOL() )
00397 {
00398
00399 AssertError[nAsserts] = ErrorDefault;
00400 }
00401
00402 lgQuantityLog[nAsserts] = true;
00403 }
00404
00405
00406 else if( p.nMatch(" LINE " ) )
00407 {
00408 if( p.nMatch("LUMI") || p.nMatch("INTE"))
00409 {
00410
00411 strcpy( chAssertType[nAsserts] , "Ll" );
00412
00413
00414 lgQuantityLog[nAsserts] = true;
00415 }
00416 else
00417 {
00418
00419 strcpy( chAssertType[nAsserts] , "Lr" );
00420
00421
00422 lgQuantityLog[nAsserts] = false;
00423 }
00424
00425
00426
00427 p.GetQuote( chLabel , true );
00428
00429
00430 if( chLabel[4] != '\0' )
00431 {
00432 fprintf( ioQQQ,
00433 " The label must be exactly 4 char long, between double quotes.\n");
00434 fprintf( ioQQQ, " Sorry.\n" );
00435 cdEXIT(EXIT_FAILURE);
00436 }
00437
00438
00439 strcpy( chAssertLineLabel[nAsserts], chLabel );
00440
00441
00442 wavelength[nAsserts] = (realnum)p.getWave();
00443
00444
00445
00446 AssertQuantity[nAsserts] = p.FFmtRead();
00447 if( p.lgEOL() )
00448 {
00449 p.NoNumb("intensity/luminosity");
00450 }
00451
00452 if( lgQuantityLog[nAsserts] )
00453 {
00454 if( AssertQuantity[nAsserts] > DBL_MAX_10_EXP ||
00455 AssertQuantity[nAsserts] < -DBL_MAX_10_EXP )
00456 {
00457 fprintf( ioQQQ,
00458 " The asserted quantity is a log, but is too large or small, value is %e.\n",
00459 AssertQuantity[nAsserts] );
00460 fprintf( ioQQQ, " I would crash if I tried to use it.\n Sorry.\n" );
00461 cdEXIT(EXIT_FAILURE);
00462 }
00463 AssertQuantity[nAsserts] =
00464 pow(10.,AssertQuantity[nAsserts] );
00465 }
00466 if( AssertQuantity[nAsserts]<= 0. )
00467 {
00468 fprintf( ioQQQ,
00469 " The relative intensity must be positive, and was %e.\n",AssertQuantity[nAsserts] );
00470 fprintf( ioQQQ, " Sorry.\n" );
00471 cdEXIT(EXIT_FAILURE);
00472 }
00473
00474
00475 AssertError[nAsserts] = p.FFmtRead();
00476 if( p.lgEOL() )
00477 {
00478
00479 AssertError[nAsserts] = ErrorDefault;
00480 }
00481 }
00482
00483
00484 else if( p.nMatch("CASE" ) )
00485 {
00486
00487 strcpy( chAssertType[nAsserts] , "CB" );
00488
00489 AssertError[nAsserts] = p.FFmtRead();
00490 if( p.lgEOL() )
00491
00492 AssertError[nAsserts] = ErrorDefault;
00493 AssertQuantity[nAsserts] = 0;
00494 wavelength[nAsserts] = 0.;
00495
00496
00497
00498 if( p.GetParam("FAINT",&Param[nAsserts][4]) )
00499 {
00500 if( p.lgEOL() )
00501 {
00502
00503 fprintf(ioQQQ," The monitor Case B faint option must have a number,"
00504 " the relative intensity of the fainest line to monitor.\n");
00505 cdEXIT(EXIT_FAILURE);
00506 }
00507
00508 if( Param[nAsserts][4]<=0. )
00509 Param[nAsserts][4] = pow(10., Param[nAsserts][4] );
00510 }
00511 else
00512 {
00513
00514 Param[nAsserts][4] = SMALLFLOAT;
00515 }
00516
00517
00518 if( p.GetRange("RANG",&Param[nAsserts][2],&Param[nAsserts][3]) )
00519 {
00520 if( p.lgEOL() )
00521 {
00522
00523 fprintf(ioQQQ," The monitor Case B range option must have two numbers,"
00524 " the lower and upper limit to the wavelengths in Angstroms.\n");
00525 fprintf(ioQQQ," There must be a total of three numbers on the line,"
00526 " the relative error followed by the lower and upper limits to the "
00527 "wavelength in Angstroms.\n");
00528 cdEXIT(EXIT_FAILURE);
00529 }
00530 if( Param[nAsserts][2]>Param[nAsserts][3])
00531 {
00532
00533 double sav = Param[nAsserts][3];
00534 Param[nAsserts][3] = Param[nAsserts][2];
00535 Param[nAsserts][2] = sav;
00536 }
00537 }
00538 else
00539 {
00540
00541 Param[nAsserts][2] = 0.;
00542 Param[nAsserts][3] = 1e30;
00543 }
00544
00545 if( p.nMatch("H-LI" ) )
00546 {
00547
00548 if( (nelem = p.GetElem()) < 0 )
00549 {
00550
00551 fprintf(ioQQQ, "monitor H-like case B did not find an element on this line, sorry\n");
00552 p.PrintLine(ioQQQ);
00553 cdEXIT(EXIT_FAILURE);
00554 }
00555 if( nelem>7 )
00556 {
00557
00558 fprintf(ioQQQ, "monitor H-like cannot do elements heavier than O, sorry\n");
00559 p.PrintLine(ioQQQ);
00560 cdEXIT(EXIT_FAILURE);
00561 }
00562 Param[nAsserts][0] = ipH_LIKE;
00563 Param[nAsserts][1] = nelem;
00564
00565 strcpy( chAssertLineLabel[nAsserts], "Ca ");
00566 if( p.nMatch("CASE A " ) )
00567 strcat( chAssertLineLabel[nAsserts] , "A");
00568 else
00569 strcat( chAssertLineLabel[nAsserts] , "B");
00570 }
00571 else if( p.nMatch("HE-L") )
00572 {
00573
00574 Param[nAsserts][0] = ipHE_LIKE;
00575 Param[nAsserts][1] = ipHELIUM;
00576
00577 strcpy( chAssertLineLabel[nAsserts] , "Ca B");
00578 }
00579 else
00580 {
00581
00582 fprintf( ioQQQ,
00583 " I could not identify an iso-sequence on this Case A/B command.\n");
00584 fprintf( ioQQQ, " Sorry.\n" );
00585 cdEXIT(EXIT_FAILURE);
00586 }
00587 }
00588
00589
00590 else if( p.nMatch("DEPA") )
00591 {
00592
00593 AssertQuantity[nAsserts] = p.FFmtRead();
00594 if( p.lgEOL() )
00595 p.NoNumb("average departure coefficient");
00596
00597
00598 AssertError[nAsserts] = p.FFmtRead();
00599 if( p.lgEOL() )
00600
00601 AssertError[nAsserts] = ErrorDefault;
00602
00603
00604 if( p.nMatch("H-LI" ) )
00605 {
00606 Param[nAsserts][0] = ipH_LIKE;
00607 strcpy( chAssertLineLabel[nAsserts] , "dHyd" );
00608
00609 strcpy( chAssertType[nAsserts] , "DI" );
00610
00611 if( (wavelength[nAsserts] = (realnum)p.GetElem()) < 0 )
00612 {
00613 fprintf( ioQQQ,
00614 " I could not identify an element name on this line.\n");
00615 fprintf( ioQQQ, " Sorry.\n" );
00616 cdEXIT(EXIT_FAILURE);
00617 }
00618 if( p.nMatch("ZEROOK") )
00619 Param[nAsserts][1] = 1.;
00620 else
00621 Param[nAsserts][1] = 0.;
00622 }
00623
00624
00625 else if( p.nMatch("HE-L" ) )
00626 {
00627 Param[nAsserts][0] = ipHE_LIKE;
00628 strcpy( chAssertLineLabel[nAsserts] , "dHel" );
00629
00630 strcpy( chAssertType[nAsserts] , "DI" );
00631
00632 if( (wavelength[nAsserts] = (realnum)p.GetElem()) < 0 )
00633 {
00634 fprintf( ioQQQ,
00635 " I could not identify an element name on this line.\n");
00636 fprintf( ioQQQ, " Sorry.\n" );
00637 cdEXIT(EXIT_FAILURE);
00638 }
00639 if( p.nMatch("ZEROOK") )
00640 Param[nAsserts][1] = 1.;
00641 else
00642 Param[nAsserts][1] = 0.;
00643 }
00644
00645
00646 else if( p.nMatch("FEII") || p.nMatch("FE II" ) )
00647 {
00648
00649 strcpy( chAssertLineLabel[nAsserts] , "d Fe" );
00650
00651 strcpy( chAssertType[nAsserts] , "d " );
00652
00653 wavelength[nAsserts] = 2;
00654 }
00655
00656
00657 else if( p.nMatch("HMIN" ) )
00658 {
00659
00660 strcpy( chAssertLineLabel[nAsserts] , "d H-" );
00661
00662 strcpy( chAssertType[nAsserts] , "d-" );
00663
00664 wavelength[nAsserts] = -1;
00665 }
00666 else
00667 {
00668 fprintf( ioQQQ,
00669 " There must be a second key: FEII, H-LIke, HMINus, or HE-Like followed by element.\n");
00670 fprintf( ioQQQ, " Sorry.\n" );
00671 cdEXIT(EXIT_FAILURE);
00672 }
00673
00674
00675 if( p.nMatch("EXCI" ) )
00676 {
00677
00678 AssertQuantity2[nAsserts] = 1.;
00679 }
00680 else
00681 {
00682
00683 AssertQuantity2[nAsserts] = 0.;
00684 }
00685 }
00686
00687
00688 else if( p.nMatch(" MAP" ) )
00689 {
00690
00691
00692
00693 if( p.nMatch("HEAT" ) )
00694 {
00695 strcpy( chAssertType[nAsserts] , "mh" );
00696 }
00697 else if( p.nMatch("COOL" ) )
00698 {
00699 strcpy( chAssertType[nAsserts] , "mc" );
00700 }
00701 else
00702 {
00703 fprintf( ioQQQ,
00704 " There must be a second key, HEATing or COOLing.\n");
00705 fprintf( ioQQQ, " Sorry.\n" );
00706 cdEXIT(EXIT_FAILURE);
00707 }
00708
00709
00710 AssertQuantity2[nAsserts] = p.FFmtRead();
00711 if( p.lgEOL() )
00712 {
00713 p.NoNumb("temperature");
00714 }
00715
00716 if( AssertQuantity2[nAsserts] <= 10. )
00717 {
00718
00719 AssertQuantity2[nAsserts] =
00720 pow(10.,AssertQuantity2[nAsserts] );
00721 }
00722
00723
00724 wavelength[nAsserts] = (realnum)AssertQuantity2[nAsserts];
00725
00726
00727 AssertQuantity[nAsserts] = p.FFmtRead();
00728 if( p.lgEOL() )
00729 {
00730 p.NoNumb("heating/cooling");
00731 }
00732
00733
00734 AssertQuantity[nAsserts] = pow(10., AssertQuantity[nAsserts]);
00735
00736
00737
00738 AssertError[nAsserts] = p.FFmtRead();
00739 if( p.lgEOL() )
00740 {
00741
00742 AssertError[nAsserts] = ErrorDefault;
00743 }
00744
00745
00746 lgQuantityLog[nAsserts] = true;
00747 }
00748
00749
00750 else if( p.nMatch("COLU" ) )
00751 {
00752
00753 strcpy( chAssertType[nAsserts] , "cd" );
00754
00755
00756 wavelength[nAsserts] = 0;
00757
00758 char chLabel[INPUT_LINE_LENGTH];
00759
00760 if ( p.GetQuote(chLabel,false) == 0 )
00761 {
00762
00763 for (long i=strlen(chLabel);i<NCHAR-1;++i)
00764 chLabel[i] = ' ';
00765 chLabel[NCHAR-1] = '\0';
00766 strcpy( chAssertLineLabel[nAsserts], chLabel );
00767 }
00768 else if( p.nMatchErase(" H2 ") )
00769 {
00770 strcpy( chAssertLineLabel[nAsserts], "H2 " );
00771 }
00772 else if( p.nMatchErase("H3+ "))
00773 {
00774 strcpy( chAssertLineLabel[nAsserts], "H3+ " );
00775 }
00776 else if( p.nMatchErase("H2+ "))
00777 {
00778 strcpy( chAssertLineLabel[nAsserts], "H2+ " );
00779 }
00780 else if( p.nMatchErase(" H- "))
00781 {
00782 strcpy( chAssertLineLabel[nAsserts], "H- " );
00783 }
00784 else if( p.nMatchErase("H2G "))
00785 {
00786 strcpy( chAssertLineLabel[nAsserts], "H2g " );
00787 }
00788 else if( p.nMatchErase("H2* "))
00789 {
00790 strcpy( chAssertLineLabel[nAsserts], "H2* " );
00791 }
00792 else if( p.nMatchErase("HEH+"))
00793 {
00794 strcpy( chAssertLineLabel[nAsserts], "HeH+" );
00795 }
00796 else if( p.nMatchErase(" O2 "))
00797 {
00798 strcpy( chAssertLineLabel[nAsserts], "O2 " );
00799 }
00800 else if( p.nMatchErase("H2O "))
00801 {
00802 strcpy( chAssertLineLabel[nAsserts], "H2O " );
00803 }
00804 else if( p.nMatchErase(" C2 "))
00805 {
00806 strcpy( chAssertLineLabel[nAsserts], "C2 " );
00807 }
00808 else if( p.nMatchErase(" C3 "))
00809 {
00810 strcpy( chAssertLineLabel[nAsserts], "C3 " );
00811 }
00812 else if( p.nMatch(" CO "))
00813 {
00814 strcpy( chAssertLineLabel[nAsserts], "CO " );
00815 }
00816 else if( p.nMatch("SIO "))
00817 {
00818 strcpy( chAssertLineLabel[nAsserts], "SiO " );
00819 }
00820 else if( p.nMatch(" OH "))
00821 {
00822 strcpy( chAssertLineLabel[nAsserts], "OH " );
00823 }
00824 else if( p.nMatch(" CN ") )
00825 {
00826 strcpy( chAssertLineLabel[nAsserts], "CN " );
00827 }
00828 else if( p.nMatch(" CH ") )
00829 {
00830 strcpy( chAssertLineLabel[nAsserts], "CH " );
00831 }
00832 else if( p.nMatch(" CH+") )
00833 {
00834 strcpy( chAssertLineLabel[nAsserts], "CH+ " );
00835 }
00836 else
00837 {
00838 fprintf( ioQQQ,
00839 " I could not identify H2, H3+, H2+, H2g, H2*, H2H+, CO, O2, SiO, OH, C2 or C3 or on this line.\n");
00840 fprintf( ioQQQ, " Sorry.\n" );
00841 cdEXIT(EXIT_FAILURE);
00842 }
00843
00844 if( p.nMatch( "LEVE" ) )
00845 {
00846 if (strncmp( chAssertLineLabel[nAsserts], "H2 ", NCHAR) != 0)
00847 {
00848 fprintf( ioQQQ, " LEVEL option only implemented for H2.\n");
00849 fprintf( ioQQQ, " Sorry.\n" );
00850 cdEXIT(EXIT_FAILURE);
00851 }
00852
00853
00854 Param[nAsserts][0] = p.FFmtRead();
00855 if( p.lgEOL() )
00856 p.NoNumb("level v" );
00857 Param[nAsserts][1] = p.FFmtRead();
00858 if( p.lgEOL() )
00859 p.NoNumb("level J" );
00860
00861 wavelength[nAsserts] = (realnum)(100.*Param[nAsserts][0] + Param[nAsserts][1]);
00862 }
00863 else
00864 {
00865
00866 Param[nAsserts][0] = -1.;
00867 Param[nAsserts][1] = -1.;
00868 }
00869
00870
00871
00872 AssertQuantity[nAsserts] = p.FFmtRead();
00873 if( p.lgEOL() )
00874 {
00875 p.NoNumb("column density");
00876 }
00877
00878 AssertQuantity[nAsserts] =
00879 pow(10.,AssertQuantity[nAsserts] );
00880
00881
00882
00883 AssertError[nAsserts] = p.FFmtRead();
00884 if( p.lgEOL() )
00885 {
00886
00887 AssertError[nAsserts] = ErrorDefault;
00888 }
00889
00890
00891
00892 if( p.nMatch( " LOG" ) )
00893 {
00894 AssertError[nAsserts] = pow(10. , AssertError[nAsserts] );
00895 }
00896
00897
00898 lgQuantityLog[nAsserts] = true;
00899 }
00900
00901
00902 else if( p.nMatch("GRAI") && p.nMatch(" H2 ") )
00903 {
00904
00905 strcpy( chAssertType[nAsserts] , "g2" );
00906
00907 strcpy( chAssertLineLabel[nAsserts], "R H2" );
00908
00909 nelem = (long int)p.FFmtRead();
00910 if( nelem!=2 )
00911 {
00912 fprintf( ioQQQ,
00913 " I did not find a 2, as in H2, as the first number on this line.\n");
00914 fprintf( ioQQQ,
00915 " The rate should be the second number.\n");
00916 fprintf( ioQQQ, " Sorry.\n" );
00917 cdEXIT(EXIT_FAILURE);
00918 }
00919
00920 AssertQuantity[nAsserts] = p.FFmtRead();
00921 if( p.lgEOL() )
00922 {
00923 p.NoNumb("grain H2 formation");
00924 }
00925
00926 if( AssertQuantity[nAsserts] <=0. )
00927 AssertQuantity[nAsserts] = pow(10.,AssertQuantity[nAsserts] );
00928
00929 wavelength[nAsserts] = 0;
00930
00931
00932
00933 AssertError[nAsserts] = p.FFmtRead();
00934 if( p.lgEOL() )
00935 {
00936
00937 AssertError[nAsserts] = ErrorDefault;
00938
00939 lgQuantityLog[nAsserts] = true;
00940 }
00941 }
00942
00943
00944 else if( p.nMatch( "GRAI" ) && p.nMatch( "POTE") )
00945 {
00946
00947 strcpy( chAssertType[nAsserts] , "gp" );
00948
00949 strcpy( chAssertLineLabel[nAsserts], "GPot" );
00950
00951
00952 wavelength[nAsserts] = (realnum)p.FFmtRead();
00953
00954 AssertQuantity[nAsserts] = p.FFmtRead();
00955
00956 if( p.lgEOL() )
00957 {
00958 p.NoNumb("grain potential");
00959 }
00960
00961
00962
00963 AssertError[nAsserts] = p.FFmtRead();
00964
00965 if( p.lgEOL() )
00966 {
00967
00968 AssertError[nAsserts] = ErrorDefault;
00969 }
00970 }
00971
00972
00973 else if( p.nMatch("TEMP") )
00974 {
00975
00976
00977
00978
00979
00980 if( p.nMatch("VOLU") )
00981 {
00982 strcpy( chAssertType[nAsserts] , "T ");
00983 }
00984 else
00985 {
00986
00987 strcpy( chAssertType[nAsserts] , "t ");
00988 }
00989
00990
00991
00992 if( p.nMatch("GRAI") )
00993 {
00994 strcpy( chAssertLineLabel[nAsserts], "GTem" );
00995
00996
00997 nelem = LONG_MAX-2;
00998
00999
01000
01001
01002
01003 }
01004
01005
01006 else if( p.nMatch("FACE") )
01007 {
01008 strcpy( chAssertLineLabel[nAsserts], "face" );
01009
01010 nelem = 0;
01011 }
01012
01013
01014 else if( p.nMatch(" H2 ") )
01015 {
01016 strcpy( chAssertLineLabel[nAsserts], "H2 " );
01017
01018 nelem = 0;
01019 }
01020
01021
01022
01023 else if( (nelem = p.GetElem()) >= 0 )
01024 {
01025
01026
01027 strcpy( chAssertLineLabel[nAsserts], elementnames.chElementNameShort[nelem] );
01028 }
01029 else
01030 {
01031 fprintf( ioQQQ,
01032 " I could not identify an element name on this line.\n");
01033 fprintf( ioQQQ, " Sorry.\n" );
01034 cdEXIT(EXIT_FAILURE);
01035 }
01036
01037
01038 if( strcmp( chAssertLineLabel[nAsserts], "face" )== 0)
01039 {
01040
01041 wavelength[nAsserts] = 0;
01042 }
01043 else if( strcmp( chAssertLineLabel[nAsserts], "H2 " )== 0)
01044 {
01045 int j;
01046
01047 wavelength[nAsserts] = 0;
01048
01049 j = (int)p.FFmtRead();
01050 if( j != 2 )
01051 TotalInsanity();
01052 ++i;
01053 }
01054 else
01055 {
01056 wavelength[nAsserts] = (realnum)p.FFmtRead();
01057 if( p.lgEOL() )
01058 {
01059 p.NoNumb("ionization stage");
01060 }
01061
01062 if( wavelength[nAsserts] < 1 || wavelength[nAsserts] > nelem+2 )
01063 {
01064 fprintf( ioQQQ,
01065 " The ionization stage is inappropriate for this element.\n");
01066 fprintf( ioQQQ, " Sorry.\n" );
01067 cdEXIT(EXIT_FAILURE);
01068 }
01069 }
01070
01071 if( p.nMatch( "GRID" ) )
01072 {
01073 AssertQuantity[nAsserts] = AssertVector[optimize.nOptimiz];
01074 }
01075 else
01076 {
01077
01078 AssertQuantity[nAsserts] = p.FFmtRead();
01079 if( p.lgEOL() )
01080 p.NoNumb("temperature");
01081 }
01082
01083
01084 AssertError[nAsserts] = p.FFmtRead();
01085 if( p.lgEOL() )
01086 AssertError[nAsserts] = ErrorDefault;
01087
01088
01089
01090 if( AssertQuantity[nAsserts] <= 10. && !p.nMatch( "LINE" ) )
01091 {
01092
01093 AssertQuantity[nAsserts] =
01094 pow(10.,AssertQuantity[nAsserts] );
01095
01096 lgQuantityLog[nAsserts] = true;
01097 }
01098 }
01099
01100
01101 else if( p.nMatch("HHEI") )
01102 {
01103
01104 strcpy( chAssertType[nAsserts] , "hh" );
01105
01106 strcpy( chAssertLineLabel[nAsserts], "HHei" );
01107
01108
01109
01110 AssertQuantity[nAsserts] = p.FFmtRead();
01111 if( p.lgEOL() )
01112 {
01113 p.NoNumb("ionization correction factor");
01114 }
01115
01116 wavelength[nAsserts] = 0;
01117
01118
01119
01120 AssertError[nAsserts] = p.FFmtRead();
01121 if( p.lgEOL() )
01122 {
01123
01124 AssertError[nAsserts] = ErrorDefault;
01125 }
01126 }
01127
01128
01129 else if( p.nMatch(" H2 ") )
01130 {
01131
01132 if( p.nMatch("ORTH") )
01133 {
01134
01135 strcpy( chAssertType[nAsserts] , "or" );
01136
01137 strcpy( chAssertLineLabel[nAsserts], "orth" );
01138 }
01139 else
01140 {
01141 fprintf( ioQQQ,
01142 " I could not identify a second keyword on this line.\n");
01143 fprintf( ioQQQ, " Sorry.\n" );
01144 cdEXIT(EXIT_FAILURE);
01145 }
01146
01147
01148 n2 = (long)p.FFmtRead();
01149 if( p.lgEOL() || n2 != 2 )
01150 {
01151 p.NoNumb("the 2 in H2 ?!");
01152 }
01153 AssertQuantity[nAsserts] = p.FFmtRead();
01154
01155 wavelength[nAsserts] = 0;
01156
01157
01158
01159 AssertError[nAsserts] = p.FFmtRead();
01160 if( p.lgEOL() )
01161 {
01162
01163 AssertError[nAsserts] = ErrorDefault;
01164 }
01165 }
01166
01167
01168 else if( p.nMatch(" MPI") )
01169 {
01170
01171 strcpy( chAssertType[nAsserts] , "mp" );
01172
01173 strcpy( chAssertLineLabel[nAsserts], "mpi " );
01174
01175 wavelength[nAsserts] = 0.;
01176 AssertQuantity[nAsserts] = 1.;
01177 AssertError[nAsserts] = ErrorDefault;
01178 }
01179
01180
01181 else if( p.nMatch("NZON") )
01182 {
01183
01184 strcpy( chAssertType[nAsserts] , "z " );
01185
01186 strcpy( chAssertLineLabel[nAsserts], "zone" );
01187
01188
01189 wavelength[nAsserts] = 0.;
01190 AssertQuantity[nAsserts] = p.FFmtRead();
01191 if( p.lgEOL() )
01192 p.NoNumb("zone number");
01193
01194
01195 AssertError[nAsserts] = p.FFmtRead();
01196 if( p.lgEOL() )
01197 AssertError[nAsserts] = ErrorDefault;
01198 }
01199
01200
01201 else if( p.nMatch("PRES") && p.nMatch("ERRO") )
01202 {
01203
01204 strcpy( chAssertType[nAsserts] , "pr" );
01205
01206 strcpy( chAssertLineLabel[nAsserts], "pres" );
01207
01208
01209
01210 wavelength[nAsserts] = 0;
01211 AssertQuantity[nAsserts] = (double)p.FFmtRead();
01212 if( p.lgEOL() )
01213 {
01214 p.NoNumb("pressure error");
01215 }
01216 else if( AssertQuantity[nAsserts] <= 0.)
01217 {
01218
01219 AssertQuantity[nAsserts] = pow(10. , AssertQuantity[nAsserts]);
01220 }
01221
01222
01223
01224 AssertError[nAsserts] = p.FFmtRead();
01225 if( p.lgEOL() )
01226 {
01227
01228 AssertError[nAsserts] = ErrorDefault;
01229 }
01230 }
01231
01232 else if( p.nMatch("PRADMAX") )
01233 {
01234
01235
01236 strcpy( chAssertType[nAsserts] , "RM" );
01237
01238 strcpy( chAssertLineLabel[nAsserts], "Prad" );
01239
01240
01241
01242 wavelength[nAsserts] = 0;
01243 AssertQuantity[nAsserts] = (double)p.FFmtRead();
01244 if( p.lgEOL() )
01245 {
01246 p.NoNumb("PRADMAX");
01247 }
01248 else if( AssertQuantity[nAsserts] <= 0.)
01249 {
01250
01251 AssertQuantity[nAsserts] = pow(10. , AssertQuantity[nAsserts]);
01252 }
01253
01254
01255
01256 AssertError[nAsserts] = p.FFmtRead();
01257 if( p.lgEOL() )
01258 {
01259
01260 AssertError[nAsserts] = ErrorDefault;
01261 }
01262 }
01263
01264
01265 else if( p.nMatch("CSUP") )
01266 {
01267
01268 strcpy( chAssertType[nAsserts] , "sc" );
01269
01270 strcpy( chAssertLineLabel[nAsserts], "sion" );
01271
01272
01273 AssertQuantity[nAsserts] = p.FFmtRead();
01274 if( p.lgEOL() )
01275 {
01276 p.NoNumb("secondary ionization rate");
01277 }
01278
01279 AssertQuantity[nAsserts] = pow(10., AssertQuantity[nAsserts] );
01280
01281
01282 wavelength[nAsserts] = 0;
01283
01284
01285
01286 AssertError[nAsserts] = p.FFmtRead();
01287 if( p.lgEOL() )
01288 {
01289
01290 AssertError[nAsserts] = ErrorDefault;
01291 }
01292
01293
01294 lgQuantityLog[nAsserts] = true;
01295 }
01296
01297
01298 else if( p.nMatch("HTOT") )
01299 {
01300
01301 strcpy( chAssertType[nAsserts] , "ht" );
01302
01303
01304 strcpy( chAssertLineLabel[nAsserts], "heat" );
01305
01306
01307 AssertQuantity[nAsserts] = p.FFmtRead();
01308 if( p.lgEOL() )
01309 {
01310 p.NoNumb("heating rate");
01311 }
01312
01313 AssertQuantity[nAsserts] = pow(10., AssertQuantity[nAsserts] );
01314
01315
01316 wavelength[nAsserts] = 0;
01317
01318
01319
01320 AssertError[nAsserts] = p.FFmtRead();
01321 if( p.lgEOL() )
01322 {
01323
01324 AssertError[nAsserts] = ErrorDefault;
01325 }
01326
01327
01328 lgQuantityLog[nAsserts] = true;
01329 }
01330
01331
01332 else if( p.nMatch("CTOT") )
01333 {
01334
01335 strcpy( chAssertType[nAsserts] , "ct" );
01336
01337
01338 strcpy( chAssertLineLabel[nAsserts], "cool" );
01339
01340
01341 if( p.nMatch( "GRID" ) )
01342 {
01343 AssertQuantity[nAsserts] = AssertVector[optimize.nOptimiz];
01344 }
01345 else
01346 {
01347 AssertQuantity[nAsserts] = p.FFmtRead();
01348
01349 if( p.lgEOL() )
01350 {
01351 p.NoNumb("cooling rate");
01352 }
01353 }
01354
01355 AssertQuantity[nAsserts] = pow(10., AssertQuantity[nAsserts] );
01356
01357
01358 wavelength[nAsserts] = 0;
01359
01360
01361
01362 AssertError[nAsserts] = p.FFmtRead();
01363 if( p.lgEOL() )
01364 {
01365
01366 AssertError[nAsserts] = ErrorDefault;
01367 }
01368
01369
01370 lgQuantityLog[nAsserts] = true;
01371 }
01372
01373
01374 else if( p.nMatch("ITRZ") )
01375 {
01376
01377 strcpy( chAssertType[nAsserts] , "iz" );
01378
01379 strcpy( chAssertLineLabel[nAsserts], "itrz" );
01380
01381
01382 AssertQuantity[nAsserts] = p.FFmtRead();
01383 if( p.lgEOL() )
01384 {
01385 p.NoNumb("iterations per zone");
01386 }
01387
01388 wavelength[nAsserts] = 0;
01389
01390
01391 AssertError[nAsserts] = p.FFmtRead();
01392 if( p.lgEOL() )
01393 AssertError[nAsserts] = ErrorDefaultPerformance;
01394 }
01395
01396
01397 else if( p.nMatch("EDEN") )
01398 {
01399
01400 strcpy( chAssertType[nAsserts] , "e " );
01401
01402 strcpy( chAssertLineLabel[nAsserts], "eden" );
01403
01404
01405 AssertQuantity[nAsserts] =
01406 pow(10., p.FFmtRead() );
01407 if( p.lgEOL() )
01408 {
01409 p.NoNumb(" electron density of the last zone");
01410 }
01411
01412
01413
01414 AssertError[nAsserts] = p.FFmtRead();
01415 if( p.lgEOL() )
01416 {
01417
01418 AssertError[nAsserts] = ErrorDefault;
01419 }
01420 wavelength[nAsserts] = 0;
01421
01422
01423 lgQuantityLog[nAsserts] = true;
01424 }
01425
01426
01427 else if( p.nMatch(" TU ") )
01428 {
01429
01430 strcpy( chAssertType[nAsserts] , "Tu" );
01431 strcpy( chAssertLineLabel[nAsserts], "Tu " );
01432
01433
01434 AssertQuantity[nAsserts] = p.FFmtRead();
01435 if( p.lgEOL() )
01436 p.NoNumb("energy density of last zone");
01437
01438
01439
01440 lgQuantityLog[nAsserts] = false;
01441 if( AssertQuantity[nAsserts] <= 10. && !p.nMatch( "LINE" ) )
01442 {
01443
01444 AssertQuantity[nAsserts] =
01445 pow(10.,AssertQuantity[nAsserts] );
01446
01447 }
01448
01449
01450
01451 AssertError[nAsserts] = p.FFmtRead();
01452 if( p.lgEOL() )
01453
01454 AssertError[nAsserts] = ErrorDefault;
01455 wavelength[nAsserts] = 0;
01456 }
01457
01458
01459 else if( p.nMatch("THIC") || p.nMatch("DEPT") )
01460 {
01461
01462 strcpy( chAssertType[nAsserts] , "th" );
01463
01464 strcpy( chAssertLineLabel[nAsserts], "thic" );
01465
01466
01467 AssertQuantity[nAsserts] = pow(10., p.FFmtRead() );
01468 if( p.lgEOL() )
01469 {
01470 p.NoNumb("thickness or depth of model");
01471 }
01472
01473
01474
01475 AssertError[nAsserts] = p.FFmtRead();
01476 if( p.lgEOL() )
01477 {
01478
01479 AssertError[nAsserts] = ErrorDefault;
01480 }
01481 wavelength[nAsserts] = 0;
01482
01483
01484 lgQuantityLog[nAsserts] = true;
01485 }
01486
01487
01488 else if( p.nMatch("RADI") )
01489 {
01490
01491 strcpy( chAssertType[nAsserts] , "ra" );
01492
01493 strcpy( chAssertLineLabel[nAsserts], "radi" );
01494
01495
01496 AssertQuantity[nAsserts] = pow(10., p.FFmtRead() );
01497 if( p.lgEOL() )
01498 {
01499 p.NoNumb("outer radius");
01500 }
01501
01502
01503
01504 AssertError[nAsserts] = p.FFmtRead();
01505 if( p.lgEOL() )
01506 {
01507
01508 AssertError[nAsserts] = ErrorDefault;
01509 }
01510 wavelength[nAsserts] = 0;
01511
01512
01513 lgQuantityLog[nAsserts] = true;
01514 }
01515
01516
01517 else if( p.nMatch("NITE") )
01518 {
01519
01520 strcpy( chAssertType[nAsserts] , "Z " );
01521
01522 strcpy( chAssertLineLabel[nAsserts], "iter" );
01523
01524 wavelength[nAsserts] = 0.f;
01525
01526
01527 AssertQuantity[nAsserts] = p.FFmtRead();
01528 if( p.lgEOL() )
01529 {
01530 p.NoNumb("number of iterations");
01531 }
01532
01533
01534
01535 AssertError[nAsserts] = p.FFmtRead();
01536 if( p.lgEOL() )
01537 {
01538
01539 AssertError[nAsserts] = ErrorDefault;
01540 }
01541 }
01542
01543
01544
01545 else if( p.nMatch("VELO") )
01546 {
01547
01548 strcpy( chAssertType[nAsserts] , "Ve" );
01549
01550 strcpy( chAssertLineLabel[nAsserts], "vel " );
01551 wavelength[nAsserts] = 0;
01552
01553
01554 AssertQuantity[nAsserts] = p.FFmtRead();
01555 if( p.lgEOL() )
01556 p.NoNumb("terminal velocity");
01557
01558
01559
01560 AssertError[nAsserts] = p.FFmtRead();
01561 if( p.lgEOL() )
01562 {
01563
01564 AssertError[nAsserts] = ErrorDefault;
01565 }
01566 }
01567
01568 else if( p.nMatch("NOTH") )
01569 {
01570 strcpy( chAssertType[nAsserts] , "NO" );
01571 strcpy( chAssertLineLabel[nAsserts], "noth" );
01572 wavelength[nAsserts] = 0;
01573 AssertQuantity[nAsserts] = 0.;
01574 AssertError[nAsserts] = ErrorDefault;
01575 }
01576 else
01577 {
01578
01579 fprintf( ioQQQ,
01580 " Unrecognized command. The line image was\n");
01581 p.PrintLine(ioQQQ);
01582 fprintf( ioQQQ,
01583 " The options I know about are: ionization, line, departure coefficient, map, column, "
01584 "temperature, nzone, csupre, htot, itrz, eden, thickness, niter, \n");
01585 fprintf( ioQQQ, " Sorry.\n" );
01586 cdEXIT(EXIT_FAILURE);
01587 }
01588
01589
01590 ++nAsserts;
01591 if( nAsserts >= NASSERTS )
01592 {
01593 fprintf(ioQQQ,
01594 " ParseMonitorResults: too many asserts, limit is NASSERT=%d\n",
01595 NASSERTS );
01596 cdEXIT(EXIT_FAILURE);
01597 }
01598 return;
01599 }
01600
01601
01602
01603 bool lgCheckMonitors(
01604
01605
01606 FILE * ioMONITOR )
01607 {
01608 double PredQuan[NASSERTS] , RelError[NASSERTS];
01609
01610
01611 bool lgFound[NASSERTS];
01612 double relint , absint;
01613 bool lg1OK[NASSERTS];
01614 long i,j;
01615
01616
01617
01618
01619
01620
01621 multi_arr<long,2> ipDisambiguate(NASSERTS,3);
01622 long lgDisambiguate = false;
01623 char chCaps[5], chFind[5];
01624 realnum errorwave;
01625
01626 DEBUG_ENTRY( "lgCheckMonitors()" );
01627
01628
01629
01630 lgMonitorsOK = true;
01631
01632
01633 lgBigBotch = false;
01634
01635
01636
01637
01638
01639
01640
01641 if( !called.lgTalk && optimize.lgOptimize )
01642 {
01643
01644 return true;
01645 }
01646
01647
01648
01649
01650
01651 if( lines_table() )
01652 {
01653 lgBigBotch = true;
01654 lgMonitorsOK = false;
01655 }
01656
01657
01658 for(i=0; i<nAsserts; ++i )
01659 {
01660 lg1OK[i] = true;
01661 PredQuan[i] = 0.;
01662 RelError[i] = 0.;
01663 for(j=0; j<3; ++j )
01664 ipDisambiguate[i][j] = -1;
01665
01666
01667 lgFound[i] = true;
01668
01669
01670
01671 if( strcmp( chAssertType[i] , "Lr")==0 )
01672 {
01673
01674 ipDisambiguate[i][0] = cdLine( chAssertLineLabel[i] , wavelength[i] , &relint,&absint );
01675 if( ipDisambiguate[i][0] <= 0 )
01676 {
01677 fprintf( ioMONITOR,
01678 " monitor error: lgCheckMonitors could not find a line with label %s ",
01679 chAssertLineLabel[i] );
01680 prt_wl( ioMONITOR, wavelength[i] );
01681 fprintf( ioMONITOR, "\n" );
01682
01683 fprintf( ioMONITOR,
01684 " monitor error: The \"save line labels\" command is a good way to get a list of line labels.\n\n");
01685
01686 lg1OK[i] = false;
01687 RelError[i] = 100000.;
01688 PredQuan[i] = 0;
01689 lg1OK[i] = false;
01690 lgMonitorsOK = false;
01691 lgFound[i] = false;
01692 continue;
01693 }
01694 else
01695 {
01696
01697
01698
01699
01700
01701 cap4(chFind, chAssertLineLabel[i]);
01702
01703
01704 errorwave = WavlenErrorGet( wavelength[i] );
01705
01706
01707 for( j=1; j < LineSave.nsum; j++ )
01708 {
01709
01710 if( j==ipDisambiguate[i][0] )
01711 continue;
01712
01713
01714 cap4(chCaps, LineSv[j].chALab);
01715
01716
01717
01718
01719
01720
01721
01722 if( fabs(LineSv[j].wavelength-wavelength[i]) < 3.f*errorwave )
01723 {
01724
01725 if( strcmp(chCaps,chFind) == 0 )
01726 {
01727 double relint1, absint1, current_error;
01728
01729 cdLine_ip( j, &relint1, &absint1 );
01730 current_error = fabs(1. - relint1/AssertQuantity[i]);
01731
01732 if( current_error < 2.*AssertError[i] ||
01733 current_error < 2.*fabs(RelError[i]) )
01734 {
01735 lgDisambiguate = true;
01736
01737
01738 if( ipDisambiguate[i][1] > 0 )
01739 {
01740 ipDisambiguate[i][2] = j;
01741 break;
01742 }
01743 else
01744 {
01745 ipDisambiguate[i][1] = j;
01746 }
01747 }
01748 }
01749 }
01750 }
01751 }
01752
01753 PredQuan[i] = relint;
01754 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
01755 }
01756
01757
01758 else if( strcmp( chAssertType[i] , "Ll")==0 )
01759 {
01760
01761 if( cdLine( chAssertLineLabel[i] , wavelength[i] , &relint,&absint )<=0 )
01762 {
01763 fprintf( ioMONITOR,
01764 " monitor error: lgCheckMonitors could not find a line with label %s ",
01765 chAssertLineLabel[i] );
01766 prt_wl( ioMONITOR, wavelength[i] );
01767 fprintf( ioMONITOR, "\n" );
01768
01769 fprintf( ioMONITOR,
01770 " monitor error: The \"save line labels\" command is a good way to get a list of line labels.\n\n");
01771
01772 lg1OK[i] = false;
01773 RelError[i] = 10000000.;
01774 PredQuan[i] = 0;
01775 lg1OK[i] = false;
01776 lgFound[i] = false;
01777 lgMonitorsOK = false;
01778 continue;
01779 }
01780
01781 PredQuan[i] = pow(10.,absint);
01782 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
01783
01784 }
01785 else if( strcmp( chAssertType[i] , "hh" ) == 0)
01786 {
01787 double hfrac , hefrac;
01788
01789 if( cdIonFrac(
01790
01791 "hydr",
01792
01793 1,
01794
01795 &hfrac,
01796
01797 "VOLUME" ,
01798
01799 false) )
01800 {
01801 fprintf( ioMONITOR,
01802 " monitor error: lgCheckMonitors could not find h ionization fraction \n");
01803 lg1OK[i] = false;
01804 RelError[i] = 0;
01805 PredQuan[i] = 0;
01806 lgFound[i] = false;
01807 continue;
01808 }
01809 if( cdIonFrac(
01810
01811 "heli",
01812
01813 1,
01814
01815 &hefrac,
01816
01817 "VOLUME" ,
01818
01819 false) )
01820 {
01821 fprintf( ioMONITOR,
01822 " monitor error: lgCheckMonitors could not find h ionization fraction \n");
01823 lg1OK[i] = false;
01824 RelError[i] = 0;
01825 PredQuan[i] = 0;
01826 lgFound[i] = false;
01827 continue;
01828 }
01829
01830 PredQuan[i] = hefrac-hfrac;
01831
01832 RelError[i] = fabs(AssertQuantity[i] - (hefrac-hfrac) );
01833 }
01834
01835 else if( strcmp( chAssertType[i] , "mp" ) == 0)
01836 {
01837 PredQuan[i] = cpu.i().lgMPI() ? 1. : 0.;
01838
01839 RelError[i] = AssertQuantity[i] - PredQuan[i];
01840 }
01841
01842 else if( strcmp( chAssertType[i] , "z " ) == 0)
01843 {
01844
01845 PredQuan[i] = (double)nzone;
01846 if( t_version::Inst().lgRelease || t_version::Inst().lgReleaseBranch )
01847 RelError[i] = ForcePass(chAssertLimit[i]);
01848 else
01849
01850 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
01851 }
01852
01853 else if( strcmp( chAssertType[i] , "or" ) == 0)
01854 {
01855
01856 PredQuan[i] = h2.ortho_density / SDIV( h2.para_density );
01857
01858
01859 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01860 }
01861
01862 else if( strcmp( chAssertType[i] , "g2" ) == 0)
01863 {
01864
01865 PredQuan[i] = gv.rate_h2_form_grains_used_total;
01866
01867 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01868 }
01869
01870 else if( strcmp( chAssertType[i] , "RM" ) == 0)
01871 {
01872
01873 PredQuan[i] = pressure.RadBetaMax;
01874
01875 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01876 }
01877
01878 else if( strcmp( chAssertType[i] , "pr" ) == 0)
01879 {
01880
01881 double sumx=0., sumx2=0., average;
01882 long int n;
01883
01884 for( n=0; n<nzone; n++ )
01885 {
01886 sumx += struc.pressure[n];
01887 sumx2 += POW2(struc.pressure[n]);
01888 }
01889 if( nzone>1 )
01890 {
01891
01892 average = sumx/nzone;
01893
01894 sumx = sqrt( (sumx2-POW2(sumx)/nzone)/(nzone-1) );
01895
01896 PredQuan[i] = sumx / average;
01897 }
01898 else
01899 {
01900 PredQuan[i] = 0.;
01901 }
01902
01903
01904 RelError[i] = PredQuan[i];
01905 }
01906
01907 else if( strcmp( chAssertType[i] , "iz") == 0 )
01908 {
01909
01910 if( nzone > 0 )
01911 PredQuan[i] = (double)(conv.nTotalIoniz-conv.nTotalIoniz_start)/(double)(nzone);
01912 else
01913
01914 PredQuan[i] = 1e10;
01915
01916 if( t_version::Inst().lgRelease || t_version::Inst().lgReleaseBranch )
01917 RelError[i] = ForcePass(chAssertLimit[i]);
01918 else
01919
01920 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01921 }
01922
01923 else if( strcmp( chAssertType[i] , "e ") == 0 )
01924 {
01925
01926 PredQuan[i] = dense.eden;
01927
01928 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01929 }
01930
01931 else if( strcmp( chAssertType[i] , "Tu") == 0 )
01932 {
01933
01934 PredQuan[i] = pow((rfield.EnergyIncidCont+rfield.EnergyDiffCont)/
01935 SPEEDLIGHT/7.56464e-15,0.25);
01936
01937 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01938 }
01939
01940 else if( strcmp( chAssertType[i] , "th") == 0 )
01941 {
01942
01943 PredQuan[i] = radius.depth;
01944
01945 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01946 }
01947
01948 else if( strcmp( chAssertType[i] , "ra") == 0 )
01949 {
01950
01951 PredQuan[i] = radius.Radius;
01952
01953 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01954 }
01955
01956 else if( strcmp( chAssertType[i] , "Ve") == 0 )
01957 {
01958
01959 PredQuan[i] = wind.windv/1e5;
01960
01961 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01962 }
01963
01964 else if( strcmp( chAssertType[i] , "NO") == 0 )
01965 {
01966
01967 PredQuan[i] = 0;
01968
01969 RelError[i] = 0.;
01970 }
01971
01972 else if( strcmp( chAssertType[i] , "sc") == 0 )
01973 {
01974
01975 PredQuan[i] = secondaries.csupra[ipHYDROGEN][0];
01976
01977 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01978 }
01979
01980 else if( strcmp( chAssertType[i] , "ht") == 0 )
01981 {
01982
01983 PredQuan[i] = thermal.htot;
01984
01985 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01986 }
01987
01988 else if( strcmp( chAssertType[i] , "ct") == 0 )
01989 {
01990
01991 PredQuan[i] = thermal.ctot;
01992
01993 RelError[i] = 1. - PredQuan[i] / AssertQuantity[i];
01994 }
01995
01996 else if( strcmp( chAssertType[i] , "Z ") == 0 )
01997 {
01998
01999 PredQuan[i] = (double)iteration;
02000 if( t_version::Inst().lgRelease || t_version::Inst().lgReleaseBranch )
02001 RelError[i] = ForcePass(chAssertLimit[i]);
02002 else
02003
02004 RelError[i] = (double)(AssertQuantity[i] - PredQuan[i]);
02005 }
02006
02007 else if( strcmp( chAssertType[i] , "CB") == 0 )
02008 {
02009 long int nISOCaseB = (long)Param[i][0];
02010 long int nelemCaseB = (long)Param[i][1];
02011 char chElemLabelCaseB[5];
02012
02013 strcpy( chElemLabelCaseB , elementnames.chElementSym[nelemCaseB] );
02014 char chNumb[4];
02015 sprintf( chNumb , "%2li" , nelemCaseB+1-nISOCaseB );
02016 strcat( chElemLabelCaseB , chNumb );
02017
02018
02019 int iCase;
02020 if( strcmp( "Ca A", chAssertLineLabel[i] ) == 0 )
02021 iCase = 0;
02022 else if( strcmp( "Ca B", chAssertLineLabel[i] ) == 0 )
02023 iCase = 1;
02024 else
02025 TotalInsanity();
02026
02027 RelError[i] = 0.;
02028 long nHighestPrinted = iso_sp[nISOCaseB][nelemCaseB].n_HighestResolved_max;
02029 if( nISOCaseB == ipH_LIKE )
02030 {
02031 fprintf(ioMONITOR," Species nHi nLo Wl Computed Asserted error\n");
02032
02033
02034
02035
02036 for( long int ipLo=1+iCase; ipLo< MIN2(10,nHighestPrinted-1); ++ipLo )
02037 {
02038 for( long int ipHi=ipLo+1; ipHi< MIN2(25,nHighestPrinted); ++ipHi )
02039 {
02040
02041 realnum wl = atmdat.WaveLengthCaseB[nelemCaseB][ipHi][ipLo];
02042
02043 if( wl < Param[i][2] || wl > Param[i][3] )
02044 continue;
02045
02046 double relint , absint,CBrelint , CBabsint;
02047
02048 cdLine( chAssertLineLabel[i] , wl , &CBrelint , &CBabsint );
02049 if( CBrelint < Param[i][4] )
02050 continue;
02051 CBabsint = pow( 10., CBabsint );
02052 double error;
02053
02054 if( cdLine( chElemLabelCaseB , wl , &relint , &absint ) >0)
02055 {
02056 absint = pow( 10., absint );
02057 error = (CBabsint - absint)/MAX2(CBabsint , absint);
02058 double RelativeError = fabs(error) / AssertError[i];
02059
02060 if( RelativeError < 1. )
02061 {
02062 if( RelativeError < 0.25 )
02063 {
02064 fprintf( ioMONITOR, " ChkMonitor ");
02065 }
02066 else if( RelativeError < 0.50 )
02067 {
02068 fprintf( ioMONITOR, " ChkMonitor - ");
02069 }
02070 else if( RelativeError < 0.75 )
02071 {
02072 fprintf( ioMONITOR, " ChkMonitor -- ");
02073 }
02074 else if( RelativeError < 0.90 )
02075 {
02076 fprintf( ioMONITOR, " ChkMonitor --- ");
02077 }
02078 else if( RelativeError < 0.95 )
02079 {
02080 fprintf( ioMONITOR, " ChkMonitor ---- ");
02081 }
02082 else if( RelativeError < 0.98 )
02083 {
02084 fprintf( ioMONITOR, " ChkMonitor ----- ");
02085 }
02086 else
02087 {
02088 fprintf( ioMONITOR, " ChkMonitor ------ ");
02089 }
02090
02091 }
02092 else
02093 {
02094 fprintf( ioMONITOR, " ChkMonitor botch>>");
02095 }
02096 fprintf(ioMONITOR," %s %3li %3li ",
02097 chElemLabelCaseB , ipHi , ipLo );
02098 prt_wl(ioMONITOR, wl );
02099 fprintf(ioMONITOR," %.2e %.2e %10.3f",
02100 absint , CBabsint , error );
02101 }
02102 else
02103 TotalInsanity();
02104 if( fabs(error) > AssertError[i] )
02105 fprintf(ioMONITOR , " botch \n");
02106 else
02107 fprintf(ioMONITOR , "\n");
02108
02109 PredQuan[i] = 0;
02110 AssertQuantity[i] = 0.;
02111 RelError[i] = MAX2( RelError[i] , fabs(error) );
02112
02113
02114 MonitorResults.SumErrorCaseMonitor += RelError[i];
02115 ++MonitorResults.nSumErrorCaseMonitor;
02116
02117 }
02118 }
02119 fprintf(ioMONITOR,"\n");
02120 }
02121 else if( nISOCaseB == ipHE_LIKE )
02122 {
02123 if( !dense.lgElmtOn[ipHELIUM] )
02124 {
02125 fprintf(ioQQQ,"PROBLEM monitor case B for a He is requested but He is not "
02126 "included.\n");
02127 fprintf(ioQQQ,"Do not turn off He if you want to monitor its spectrum.\n");
02128 cdEXIT(EXIT_FAILURE);
02129 }
02130
02131
02132 fprintf(ioMONITOR," Wl Computed Asserted error\n");
02133 for( unsigned int ipLine=0; ipLine< atmdat.CaseBWlHeI.size(); ++ipLine )
02134 {
02135
02136 realnum wl = atmdat.CaseBWlHeI[ipLine];
02137
02138 if( wl < Param[i][2] || wl > Param[i][3] )
02139 continue;
02140 double relint , absint,CBrelint , CBabsint;
02141 cdLine( chAssertLineLabel[i] , wl , &CBrelint , &CBabsint );
02142 if( CBrelint < Param[i][4] )
02143 continue;
02144 CBabsint = pow( 10., CBabsint );
02145 double error;
02146 if( cdLine( chElemLabelCaseB , wl , &relint , &absint ) >0)
02147 {
02148 absint = pow( 10., absint );
02149 error = (CBabsint - absint)/MAX2(CBabsint , absint);
02150 double RelativeError = fabs(error) / AssertError[i];
02151
02152 if( RelativeError < 1. )
02153 {
02154 if( RelativeError < 0.25 )
02155 {
02156 fprintf( ioMONITOR, " ChkMonitor ");
02157 }
02158 else if( RelativeError < 0.50 )
02159 {
02160 fprintf( ioMONITOR, " ChkMonitor - ");
02161 }
02162 else if( RelativeError < 0.75 )
02163 {
02164 fprintf( ioMONITOR, " ChkMonitor -- ");
02165 }
02166 else if( RelativeError < 0.90 )
02167 {
02168 fprintf( ioMONITOR, " ChkMonitor --- ");
02169 }
02170 else if( RelativeError < 0.95 )
02171 {
02172 fprintf( ioMONITOR, " ChkMonitor ---- ");
02173 }
02174 else if( RelativeError < 0.98 )
02175 {
02176 fprintf( ioMONITOR, " ChkMonitor ----- ");
02177 }
02178 else
02179 {
02180 fprintf( ioMONITOR, " ChkMonitor ------ ");
02181 }
02182
02183 }
02184 else
02185 {
02186 fprintf( ioMONITOR, " ChkMonitor botch>>");
02187 }
02188 prt_wl(ioMONITOR, wl );
02189 fprintf(ioMONITOR," %.2e %.2e %10.3f",
02190 absint , CBabsint , error );
02191 }
02192 else
02193 TotalInsanity();
02194 if( fabs(error) > AssertError[i] )
02195 fprintf(ioMONITOR , " botch \n");
02196 else
02197 fprintf(ioMONITOR , "\n");
02198
02199 PredQuan[i] = 0;
02200 AssertQuantity[i] = 0.;
02201 RelError[i] = MAX2( RelError[i] , fabs(error) );
02202
02203
02204 MonitorResults.SumErrorCaseMonitor += RelError[i];
02205 ++MonitorResults.nSumErrorCaseMonitor;
02206 }
02207 fprintf(ioMONITOR,"\n");
02208 }
02209 else
02210 TotalInsanity();
02211 }
02212
02213
02214 else if( strcmp( chAssertType[i] , "DI") == 0 )
02215 {
02216
02217
02218 long ipISO = (long)Param[i][0];
02219 long nelem = (long)wavelength[i];
02220 if( !dense.lgElmtOn[nelem] )
02221 {
02222 fprintf(ioQQQ,"PROBLEM asserted element %ld is not turned on!\n",nelem);
02223 PredQuan[i] = 0.;
02224 RelError[i] = 0.;
02225 }
02226 else
02227 {
02228 RelError[i] = 0.;
02229 PredQuan[i] = 0.;
02230 long numPrintLevels = iso_sp[ipISO][nelem].numLevels_local - (long)AssertQuantity2[i];
02231 for( long n=(long)AssertQuantity2[i]; n<numPrintLevels+(long)AssertQuantity2[i]; ++n )
02232 {
02233 double relerror;
02234 relerror = 0.;
02235 PredQuan[i] += iso_sp[ipISO][nelem].fb[n].DepartCoef;
02236
02237 relerror = iso_sp[ipISO][nelem].fb[n].DepartCoef -1.;
02238 relerror = MAX2( relerror , 1. - iso_sp[ipISO][nelem].fb[n].DepartCoef );
02239 RelError[i] = MAX2( RelError[i] , relerror / AssertQuantity[i] );
02240 }
02241 PredQuan[i] /= (double)(numPrintLevels);
02242 if( fp_equal( Param[i][1], 1. ) && PredQuan[i]==0. )
02243 {
02244
02245 ASSERT( dense.xIonDense[nelem][nelem-ipISO]==0. || dense.xIonDense[nelem][nelem+1-ipISO]==0. );
02246 PredQuan[i] = AssertQuantity[i];
02247 RelError[i] = 0.;
02248 }
02249 }
02250 }
02251
02252
02253 else if( strcmp( chAssertType[i] , "d ") == 0 )
02254 {
02255 double BigError , StdDev;
02256
02257 AssertFeIIDep( &PredQuan[i] , &BigError , &StdDev );
02258
02259
02260
02261 RelError[i] = StdDev;
02262 }
02263
02264
02265 else if( strcmp( chAssertType[i] , "d-") == 0 )
02266 {
02267 PredQuan[i] = hmi.hmidep;
02268 RelError[i] = fabs( hmi.hmidep - 1.);
02269 }
02270
02271
02272 else if( (strcmp( chAssertType[i] , "f ") == 0) ||
02273 (strcmp( chAssertType[i] , "F ") == 0) )
02274 {
02275 char chWeight[7];
02276 if( strcmp( chAssertType[i] , "F ") == 0 )
02277 {
02278 strcpy( chWeight , "VOLUME" );
02279 }
02280 else
02281 {
02282
02283 strcpy( chWeight , "RADIUS" );
02284 }
02285
02286 if( cdIonFrac(
02287
02288 chAssertLineLabel[i],
02289
02290 (long)wavelength[i],
02291
02292 &relint,
02293
02294 chWeight ,
02295
02296 false) )
02297 {
02298 fprintf( ioMONITOR,
02299 " monitor error: lgCheckMonitors could not find a line with label %s %f \n",
02300 chAssertLineLabel[i] , wavelength[i] );
02301
02302 lg1OK[i] = false;
02303 RelError[i] = 0;
02304 PredQuan[i] = 0;
02305 lgFound[i] = false;
02306 continue;
02307 }
02308
02309 PredQuan[i] = relint;
02310 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
02311 }
02312
02313
02314 else if( strcmp( chAssertType[i] , "cd") == 0)
02315 {
02316
02317 if( (strcmp( chAssertLineLabel[i], "H2 " ) == 0) && (Param[i][0] >= 0.) )
02318 {
02319
02320
02321 if( (relint = cdH2_colden( (long)Param[i][0] , (long)Param[i][1] ) ) < 0. )
02322 {
02323 fprintf(ioQQQ," PROBLEM lgCheckMonitors did not find v=%li, J=%li for H2 column density.\n",
02324 (long)Param[i][0] , (long)Param[i][1] );
02325 lg1OK[i] = false;
02326 RelError[i] = 0;
02327 PredQuan[i] = 0;
02328 lgFound[i] = false;
02329 continue;
02330 }
02331 }
02332 else
02333 {
02334
02335 if( cdColm(
02336
02337 chAssertLineLabel[i],
02338
02339
02340 (long)wavelength[i],
02341
02342 &relint) )
02343 {
02344 fprintf( ioMONITOR,
02345 " monitor error: lgCheckMonitors could not find a molecule with label %s %f \n",
02346 chAssertLineLabel[i] , wavelength[i] );
02347
02348 lg1OK[i] = false;
02349 RelError[i] = 0;
02350 PredQuan[i] = 0;
02351 lgFound[i] = false;
02352 continue;
02353 }
02354 }
02355
02356 PredQuan[i] = relint;
02357 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
02358 }
02359
02360
02361 else if( (strcmp( chAssertType[i] , "MF") == 0) || (strcmp( chAssertType[i] , "mf") == 0) )
02362 {
02363
02364 relint = 0.;
02365 if( strcmp( chAssertLineLabel[i], "H2 ") ==0)
02366 {
02367
02368 if( cdColm("H2 " , 0,
02369
02370 &relint) )
02371 TotalInsanity();
02372
02373 relint = relint / colden.colden[ipCOL_HTOT];
02374 }
02375 else
02376 {
02377 fprintf( ioMONITOR,
02378 " monitor error: lgCheckMonitors could not find a molecule with label %s %f \n",
02379 chAssertLineLabel[i] , wavelength[i] );
02380
02381 lg1OK[i] = false;
02382 RelError[i] = 0;
02383 PredQuan[i] = 0;
02384 continue;
02385 }
02386
02387 PredQuan[i] = relint;
02388 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
02389 }
02390
02391
02392 else if( strcmp( chAssertType[i] , "mh") == 0 || strcmp( chAssertType[i] , "mc") == 0)
02393 {
02394
02395
02396 if( hcmap.nMapAlloc == 0 )
02397 {
02398
02399 fprintf( ioMONITOR,
02400 " monitor error: lgCheckMonitors cannot check map since map not done.\n");
02401
02402 lg1OK[i] = false;
02403 RelError[i] = 0;
02404 PredQuan[i] = 0;
02405 continue;
02406 }
02407
02408 if( AssertQuantity2[i]< hcmap.temap[0] || AssertQuantity2[i]> hcmap.temap[hcmap.nmap-1] )
02409 {
02410 fprintf( ioMONITOR,
02411 " monitor error: lgCheckMonitors cannot check map since temperature not within range.\n");
02412
02413 lg1OK[i] = false;
02414 RelError[i] = 0;
02415 PredQuan[i] = 0;
02416 continue;
02417 }
02418
02419
02420 j = 0;
02421 while( AssertQuantity2[i]>hcmap.temap[j]*1.001 && j < hcmap.nmap )
02422 {
02423 ++j;
02424 }
02425
02426
02427
02428 if( strcmp( chAssertType[i] , "mh") == 0 )
02429 {
02430
02431 PredQuan[i] = hcmap.hmap[j];
02432 strcpy(chAssertLineLabel[i],"MapH" );
02433 }
02434 else if( strcmp( chAssertType[i] , "mc") == 0)
02435 {
02436
02437 PredQuan[i] = hcmap.cmap[j];
02438 strcpy(chAssertLineLabel[i],"MapC" );
02439 }
02440 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
02441 }
02442
02443
02444 else if( (strcmp( chAssertType[i] , "t ") == 0) ||
02445 (strcmp( chAssertType[i] , "T ") == 0) )
02446 {
02447 char chWeight[7];
02448 if( strcmp( chAssertType[i] , "T ") == 0 )
02449 {
02450 strcpy( chWeight , "VOLUME" );
02451 }
02452 else
02453 {
02454
02455 strcpy( chWeight , "RADIUS" );
02456 }
02457
02458
02459 if( strcmp( chAssertLineLabel[i], "GTem" ) == 0 )
02460 {
02461 size_t nd;
02462
02463
02464 nd = (size_t)wavelength[i]-1;
02465 if( nd >= gv.bin.size() )
02466 {
02467 fprintf( ioQQQ, "Illegal grain number found: %f\n" , wavelength[i] );
02468 fprintf( ioQQQ, "Use 1 for first grain that is turned on, " );
02469 fprintf( ioQQQ, "2 for second, etc....\n" );
02470 fprintf( ioQQQ, "Old style grain numbers are not valid anymore !!\n" );
02471 cdEXIT(EXIT_FAILURE);
02472 }
02473 relint = gv.bin[nd]->avdust/radius.depth_x_fillfac;
02474 }
02475 else if( strcmp( chAssertLineLabel[i], "face" ) == 0 )
02476 {
02477
02478 relint = struc.testr[0];
02479 }
02480 else
02481 {
02482
02483 if( cdTemp(
02484
02485 chAssertLineLabel[i],
02486
02487 (long)wavelength[i],
02488
02489 &relint,
02490
02491 chWeight ) )
02492 {
02493 fprintf( ioMONITOR,
02494 " monitor error: lgCheckMonitors could not find an ion with label %s ion %li \n",
02495 chAssertLineLabel[i] , (long)wavelength[i] );
02496
02497 lg1OK[i] = false;
02498 RelError[i] = 0;
02499 PredQuan[i] = 0;
02500 lgFound[i] = false;
02501 continue;
02502 }
02503 }
02504
02505 PredQuan[i] = relint;
02506 RelError[i] = 1. - PredQuan[i]/AssertQuantity[i];
02507 }
02508
02509
02510 else if( strcmp( chAssertType[i], "gp" ) == 0 )
02511 {
02512
02513
02514 size_t nd = (size_t)wavelength[i]-1;
02515 if( nd >= gv.bin.size() )
02516 {
02517 fprintf( ioQQQ, "Illegal grain number found: %g\n" , wavelength[i] );
02518 fprintf( ioQQQ, "Use 1 for first grain that is turned on, " );
02519 fprintf( ioQQQ, "2 for second, etc....\n" );
02520 fprintf( ioQQQ, "Old style grain numbers are not valid anymore !!\n" );
02521 cdEXIT(EXIT_FAILURE);
02522 }
02523
02524
02525 PredQuan[i] = gv.bin[nd]->avdpot/radius.depth_x_fillfac;
02526
02527 RelError[i] = AssertQuantity[i] - PredQuan[i];
02528 }
02529
02530 else
02531 {
02532 fprintf( ioMONITOR,
02533 " monitor error: lgCheckMonitors received an insane chAssertType=%s, impossible\n",
02534 chAssertType[i] );
02535 ShowMe();
02536 cdEXIT(EXIT_FAILURE);
02537 }
02538
02539 if( chAssertLimit[i] == '=' )
02540 {
02541
02542 if( fabs(RelError[i]) > AssertError[i] )
02543 {
02544 lg1OK[i] = false;
02545 lgMonitorsOK = false;
02546 }
02547 }
02548 else if( chAssertLimit[i] == '<' )
02549 {
02550
02551
02552
02553
02554
02555 if( RelError[i] <= 0. )
02556 {
02557 lg1OK[i] = false;
02558 lgMonitorsOK = false;
02559 }
02560 }
02561 else if( chAssertLimit[i] == '>' )
02562 {
02563
02564
02565 if( RelError[i] >= 0. )
02566 {
02567 lg1OK[i] = false;
02568 lgMonitorsOK = false;
02569 }
02570 }
02571 }
02572
02573
02574 if( called.lgTalk && nAsserts>0 )
02575 {
02576 time_t now;
02577
02578
02579 if( lgDisambiguate )
02580 {
02581
02582 long sigfigsav = LineSave.sig_figs;
02583 double relint1, relint2, absint1;
02584
02585 LineSave.sig_figs = 6;
02586
02587 fprintf( ioMONITOR, "=============Line Disambiguation=======================================================\n" );
02588 fprintf( ioMONITOR, " Wavelengths || Intensities \n" );
02589 fprintf( ioMONITOR, "Label line match1 match2 match3 || asserted match1 match2 match3\n" );
02590
02591 for( i=0; i<nAsserts; ++i )
02592 {
02593 if( ipDisambiguate[i][1] > 0 )
02594 {
02595 fprintf( ioMONITOR , "%4s ", chAssertLineLabel[i] );
02596 prt_wl( ioMONITOR , wavelength[i] );
02597 fprintf( ioMONITOR , " " );
02598 prt_wl( ioMONITOR , LineSv[ipDisambiguate[i][0]].wavelength );
02599 fprintf( ioMONITOR , " " );
02600 prt_wl( ioMONITOR , LineSv[ipDisambiguate[i][1]].wavelength );
02601 fprintf( ioMONITOR , " " );
02602 if( ipDisambiguate[i][2] > 0 )
02603 {
02604 prt_wl( ioMONITOR , LineSv[ipDisambiguate[i][2]].wavelength );
02605 cdLine_ip( ipDisambiguate[i][2], &relint2, &absint1 );
02606 }
02607 else
02608 {
02609 fprintf( ioMONITOR , "--------" );
02610 relint2 = 0.0;
02611 }
02612 fprintf( ioMONITOR , " ||" );
02613
02614 cdLine_ip( ipDisambiguate[i][1], &relint1, &absint1 );
02615
02616 if( lgPrtSciNot )
02617 {
02618 fprintf( ioMONITOR , " %10.3e %10.3e %10.3e %10.3e\n",
02619 AssertQuantity[i],
02620 PredQuan[i] ,
02621 relint1,
02622 relint2 );
02623 }
02624 else
02625 {
02626 fprintf( ioMONITOR , " %10.4f %10.4f %10.4f %10.4f\n",
02627 AssertQuantity[i],
02628 PredQuan[i] ,
02629 relint1,
02630 relint2 );
02631 }
02632 }
02633 }
02634 fprintf( ioMONITOR, "\n" );
02635
02636
02637 LineSave.sig_figs = sigfigsav;
02638 }
02639
02640
02641 fprintf( ioMONITOR, "=============Results of monitors: Cloudy %s ", t_version::Inst().chVersion );
02642
02643
02644
02645 if( prt.lgPrintTime )
02646 {
02647
02648 now = time(NULL);
02649
02650 fprintf(ioMONITOR,"%s", ctime(&now) );
02651 }
02652 else
02653 {
02654 fprintf(ioMONITOR,"\n" );
02655 }
02656
02657 if( lgMonitorsOK )
02658 {
02659 fprintf( ioMONITOR, " No errors were found. Summary follows.\n");
02660 }
02661 else
02662 {
02663 fprintf( ioMONITOR, " Errors were found. Summary follows.\n");
02664 }
02665
02666 fprintf( ioMONITOR,
02667 " Label line computed asserted Rel Err Set err\n");
02668
02669 for( i=0; i<nAsserts; ++i )
02670 {
02671 double prtPredQuan, prtAssertQuantity;
02672
02673
02674 if( lgQuantityLog[i] )
02675 {
02676 prtPredQuan = log10( MAX2( SMALLDOUBLE , PredQuan[i] ) );
02677 prtAssertQuantity = log10( MAX2( SMALLDOUBLE , AssertQuantity[i] ) );
02678 }
02679 else
02680 {
02681 prtPredQuan = PredQuan[i];
02682 prtAssertQuantity = AssertQuantity[i];
02683 }
02684
02685 if( lg1OK[i] )
02686 {
02687
02688
02689 double relative = fabs(RelError[i] / SDIV( AssertError[i]));
02690
02691 if( relative < 0.25 || chAssertLimit[i] != '=' )
02692 {
02693 fprintf( ioMONITOR, " ChkMonitor ");
02694 }
02695 else if( relative < 0.50 )
02696 {
02697 fprintf( ioMONITOR, " ChkMonitor - ");
02698 }
02699 else if( relative < 0.75 )
02700 {
02701 fprintf( ioMONITOR, " ChkMonitor -- ");
02702 }
02703 else if( relative < 0.90 )
02704 {
02705 fprintf( ioMONITOR, " ChkMonitor --- ");
02706 }
02707 else if( relative < 0.95 )
02708 {
02709 fprintf( ioMONITOR, " ChkMonitor ---- ");
02710 }
02711 else if( relative < 0.98 )
02712 {
02713 fprintf( ioMONITOR, " ChkMonitor ----- ");
02714 }
02715 else
02716 {
02717 fprintf( ioMONITOR, " ChkMonitor ------ ");
02718 }
02719
02720 }
02721 else
02722 {
02723 fprintf( ioMONITOR, " ChkMonitor botch>>");
02724 }
02725 if( strcmp( chAssertType[i] , "Ll")==0 || ( strcmp( chAssertType[i] , "Lr")==0 ) )
02726 {
02727
02728 fprintf( ioMONITOR , "%4s ",
02729 chAssertLineLabel[i] );
02730
02731 prt_wl( ioMONITOR , wavelength[i] );
02732
02733 if( lgPrtSciNot )
02734 {
02735 fprintf( ioMONITOR , " %10.3e %c %10.3e %7.3f %7.3f ",
02736 prtPredQuan ,
02737 chAssertLimit[i] ,
02738 prtAssertQuantity ,
02739 RelError[i] ,
02740 AssertError[i]);
02741 }
02742 else
02743 {
02744 fprintf( ioMONITOR , " %10.4f %c %10.4f %7.3f %7.3f ",
02745 prtPredQuan ,
02746 chAssertLimit[i] ,
02747 prtAssertQuantity ,
02748 RelError[i] ,
02749 AssertError[i]);
02750 }
02751
02752 {
02753 enum {DEBUG_LOC=false};
02754
02755 if( DEBUG_LOC )
02756 {
02757 long ipHi, ipLo;
02758 errorwave = WavlenErrorGet( wavelength[i] );
02759
02760 for( ipLo = ipHe1s1S; ipLo <= iso_sp[ipHE_LIKE][ipHELIUM].numLevels_local -
02761 iso_sp[ipHE_LIKE][ipHELIUM].nCollapsed_local; ipLo ++ )
02762 {
02763 for( ipHi = ipLo+1; ipHi < iso_sp[ipHE_LIKE][ipHELIUM].numLevels_local -
02764 iso_sp[ipHE_LIKE][ipHELIUM].nCollapsed_local; ipHi++ )
02765 {
02766 if( fabs(iso_sp[ipHE_LIKE][ipHELIUM].trans(ipHi,ipLo).WLAng() - wavelength[i])
02767 < errorwave && ipLo!=ipHe2p3P0 && ipLo!=ipHe2p3P1 )
02768 {
02769 fprintf( ioQQQ, "\t%li %li %li %li %li %li",
02770 iso_sp[ipHE_LIKE][ipHELIUM].st[ipHi].n(),
02771 iso_sp[ipHE_LIKE][ipHELIUM].st[ipHi].l(),
02772 iso_sp[ipHE_LIKE][ipHELIUM].st[ipHi].S(),
02773 iso_sp[ipHE_LIKE][ipHELIUM].st[ipLo].n(),
02774 iso_sp[ipHE_LIKE][ipHELIUM].st[ipLo].l(),
02775 iso_sp[ipHE_LIKE][ipHELIUM].st[ipLo].S() );
02776 break;
02777 }
02778 }
02779 }
02780 }
02781 }
02782
02783 }
02784 else
02785 {
02786 fprintf( ioMONITOR , "%4s %6i %10.4f %c %10.4f %7.3f %7.3f ",
02787 chAssertLineLabel[i] ,
02788 (int)wavelength[i] ,
02789 prtPredQuan ,
02790 chAssertLimit[i] ,
02791 prtAssertQuantity ,
02792 RelError[i] ,
02793 AssertError[i]);
02794 }
02795
02796
02797
02798
02799
02800
02801 if( !lg1OK[i] && (fabs(RelError[i]) > 3.*AssertError[i]) && lgFound[i] )
02802 {
02803 fprintf( ioMONITOR , " <<BIG BOTCH!!\n");
02804 lgBigBotch = true;
02805 }
02806 else
02807 {
02808 fprintf( ioMONITOR , "\n");
02809 }
02810 }
02811 fprintf( ioMONITOR , " \n");
02812
02813
02814
02815 if( !lgMonitorsOK && lgBigBotch )
02816 {
02817
02818 fprintf( ioMONITOR, " BIG BOTCHED MONITORS!!! Big Botched Monitors!!! \n");
02819 }
02820 else if( !lgMonitorsOK )
02821 {
02822 fprintf( ioMONITOR, " BOTCHED MONITORS!!! Botched Monitors!!! \n");
02823 }
02824
02825 if( MonitorResults.nSumErrorCaseMonitor>0 )
02826 {
02827 fprintf(ioMONITOR,"\n The mean of the %li monitor Case A, B relative "
02828 "residuals is %.2f\n\n" ,
02829 MonitorResults.nSumErrorCaseMonitor,
02830 MonitorResults.SumErrorCaseMonitor /MonitorResults.nSumErrorCaseMonitor );
02831 }
02832
02833
02834 if( prt.lgPrintTime )
02835 {
02836 fprintf( ioQQQ, " %s\n\n", t_version::Inst().chInfo );
02837 }
02838 }
02839 return lgMonitorsOK;
02840 }