00001
00002
00003
00004
00005
00006
00007 #include "cddefines.h"
00008 #include "trace.h"
00009 #include "optimize.h"
00010 #include "input.h"
00011 #include "prt.h"
00012 #include "energy.h"
00013 #include "predcont.h"
00014 #include "parser.h"
00015 #include "lines_service.h"
00016
00017 static const realnum DEFERR = 0.05f;
00018
00019
00020 STATIC void GetOptLineInt(Parser &p );
00021
00022
00023 STATIC void GetOptColDen(Parser &p );
00024
00025
00026 STATIC void GetOptTemp(Parser &p );
00027
00028
00029 void ParseOptimize(
00030
00031 Parser &p)
00032 {
00033 DEBUG_ENTRY( "ParseOptimize()" );
00034
00035
00036 if( p.nMatch("FILE") )
00037 {
00038
00039
00040
00041
00042
00043
00044
00045 (void)p.GetQuote( chOptimFileName , true );
00046 }
00047
00048 else if( p.nMatch("COLU") )
00049 {
00050
00051
00052 GetOptColDen(p);
00053
00054 optimize.lgOptimize = true;
00055 }
00056
00057 else if( p.nMatch("CONT") && p.nMatch("FLUX") )
00058 {
00059
00060 double energy = p.FFmtRead();
00061 if( p.lgEOL() )
00062 p.NoNumb("energy");
00063 const char* unit = p.StandardEnergyUnit();
00064 long ind = t_PredCont::Inst().add( energy, unit );
00065 Energy E( energy, unit );
00066
00067 double flux = p.FFmtRead();
00068 if( p.lgEOL() )
00069 p.NoNumb("flux");
00070 if( flux <= 0. || p.nMatch( " LOG" ) )
00071 flux = pow(10.,flux);
00072 Flux F( E, flux, p.StandardFluxUnit() );
00073
00074 double relerr = p.FFmtRead();
00075 if( p.lgEOL() )
00076 relerr = DEFERR;
00077
00078 if( p.nMatch("<" ) )
00079 relerr = -relerr;
00080
00081 optimize.ContIndex.push_back( ind );
00082 optimize.ContEner.push_back( E );
00083 optimize.ContNFnu.push_back( F );
00084 optimize.ContNFnuErr.push_back( chi2_type(relerr) );
00085 optimize.lgOptimize = true;
00086 }
00087
00088 else if( p.nMatch("CONT") )
00089 {
00090
00091 optimize.lgOptCont = true;
00092 optimize.lgOptimize = true;
00093 }
00094
00095 else if( p.nMatch("DIAM") )
00096 {
00097
00098 optimize.optDiam = chi2_type( p.FFmtRead() );
00099 optimize.optDiamErr = chi2_type( p.FFmtRead() );
00100 if( p.lgEOL() )
00101 optimize.optDiamErr = chi2_type( DEFERR );
00102
00103 if( p.nMatch( "<" ) )
00104 optimize.optDiamErr = -optimize.optDiamErr;
00105
00106 if( p.nMatch( " LOG" ) )
00107 optimize.optDiam = chi2_type( pow(chi2_type(10.),optimize.optDiam) );
00108
00109 optimize.lgOptDiam = true;
00110
00111 optimize.lgDiamInCM = ( p.nMatch( " CM " ) != 0 );
00112 optimize.lgOptimize = true;
00113 }
00114
00115 else if( p.nMatch("INCR") )
00116 {
00117
00118 if( optimize.nparm > 0 )
00119 {
00120
00121 optimize.OptIncrm[optimize.nparm-1] =
00122 (realnum)p.FFmtRead();
00123 }
00124 }
00125
00126 else if( p.nMatch("LUMI") || p.nMatch("INTE") )
00127 {
00128
00129 optimize.optint = (realnum)p.FFmtRead();
00130 optimize.optier = (realnum)p.FFmtRead();
00131 if( p.lgEOL() )
00132 optimize.optier = DEFERR;
00133
00134
00135 optimize.nOptLum = p.nMatch("EMER") ? 1 : 0;
00136
00137
00138 optimize.lgOptLum = true;
00139 optimize.lgOptimize = true;
00140 }
00141
00142 else if( p.nMatch("ITER") )
00143 {
00144
00145 optimize.nIterOptim = (long)p.FFmtRead();
00146 }
00147
00148 else if( p.nMatch("LINE") )
00149 {
00150
00151 GetOptLineInt(p);
00152
00153 optimize.lgOptimize = true;
00154 }
00155
00156 else if( p.nMatch("PHYM") )
00157 {
00158
00159 strcpy( optimize.chOptRtn, "PHYM" );
00160 # ifdef __unix
00161
00162
00163 optimize.lgParallel = ! ( p.nMatch("SEQU") || cpu.lgMPISingleRankMode() );
00164 # else
00165 optimize.lgParallel = false;
00166 # endif
00167 if( optimize.lgParallel )
00168 {
00169 long dum = (long)p.FFmtRead();
00170
00171 optimize.useCPU = p.lgEOL() ? cpu.nCPU() : dum;
00172 }
00173 else
00174 {
00175 optimize.useCPU = 1;
00176 }
00177 }
00178
00179 else if( p.nMatch("RANG") )
00180 {
00181
00182 if( optimize.nparm > 0 )
00183 {
00184 bool lgFirstOneReal = false;
00185
00186 optimize.varang[optimize.nparm-1][0] =
00187 (realnum)p.FFmtRead();
00188 if( p.lgEOL() )
00189 {
00190 optimize.varang[optimize.nparm-1][0] = -FLT_MAX;
00191 }
00192 else
00193 {
00194 lgFirstOneReal = true;
00195 }
00196
00197 optimize.varang[optimize.nparm-1][1] =
00198 (realnum)p.FFmtRead();
00199 if( p.lgEOL() )
00200 {
00201 optimize.varang[optimize.nparm-1][1] = FLT_MAX;
00202 }
00203 else if( lgFirstOneReal )
00204 {
00205
00206
00207 ++optimize.nRangeSet;
00208 if( optimize.varang[optimize.nparm-1][1] < optimize.varang[optimize.nparm-1][0] )
00209 {
00210 realnum temp = optimize.varang[optimize.nparm-1][0];
00211 optimize.varang[optimize.nparm-1][0] = optimize.varang[optimize.nparm-1][1];
00212 optimize.varang[optimize.nparm-1][1] = temp;
00213 }
00214 }
00215 }
00216 }
00217
00218 else if( p.nMatch("SUBP") )
00219 {
00220
00221 strcpy( optimize.chOptRtn, "SUBP" );
00222 }
00223
00224
00225 else if( p.nMatch("TEMP") )
00226 {
00227
00228 GetOptTemp(p);
00229
00230 optimize.lgOptimize = true;
00231 }
00232
00233 else if( p.nMatch("TOLE") )
00234 {
00235
00236
00237 optimize.OptGlobalErr = (realnum)p.FFmtRead();
00238 }
00239
00240 else if( p.nMatch("TRAC") )
00241 {
00242 if( p.nMatch("STAR") )
00243 {
00244
00245
00246 optimize.nTrOpt = (long)p.FFmtRead();
00247 if( p.lgEOL() )
00248 {
00249 fprintf( ioQQQ, " optimize trace start command:\n" );
00250 fprintf( ioQQQ, " The iteration number must appear.\n" );
00251 cdEXIT(EXIT_FAILURE);
00252 }
00253 optimize.lgTrOpt = true;
00254 }
00255 else if( p.nMatch("FLOW") )
00256 {
00257
00258
00259 optimize.lgOptimFlow = true;
00260 }
00261 else
00262 {
00263 fprintf( ioQQQ, " optimize trace flow command:\n" );
00264 fprintf( ioQQQ, " One of the sub keys START or FLOW must appear.\n" );
00265 cdEXIT(EXIT_FAILURE);
00266 }
00267 }
00268
00269 else
00270 {
00271 p.PrintLine(ioQQQ);
00272 fprintf( ioQQQ, " is unrecognized keyword, consult HAZY.\n" );
00273 cdEXIT(EXIT_FAILURE);
00274 }
00275 return;
00276 }
00277
00278
00279 STATIC void GetOptColDen(Parser &p )
00280 {
00281
00282 DEBUG_ENTRY( "GetOptColDen()" );
00283
00284
00285
00286
00287 p.getline();
00288 if( p.m_lgEOF )
00289 {
00290 fprintf( ioQQQ, " Hit EOF while reading column density list; use END to end list.\n" );
00291 cdEXIT(EXIT_FAILURE);
00292 }
00293
00294 while( !p.m_lgEOF )
00295 {
00296
00297
00298 optimize.chColDen_label.push_back(p.getCommand(4));
00299
00300
00301
00302 long ion = nint(p.FFmtRead());
00303 if( p.lgEOL() )
00304 {
00305 p.PrintLine( ioQQQ );
00306 fprintf( ioQQQ, " The ionization stage MUST appear on this line. Sorry.\n" );
00307 cdEXIT(EXIT_FAILURE);
00308 }
00309
00310
00311
00312
00313
00314 if( ion < 0 )
00315 {
00316 p.PrintLine( ioQQQ );
00317 fprintf( ioQQQ, " An ionization stage of %ld does not make sense. Sorry.\n", ion );
00318 cdEXIT(EXIT_FAILURE);
00319 }
00320 optimize.ion_ColDen.push_back(ion);
00321
00322 optimize.ColDen_Obs.push_back( realnum(pow(10.,p.FFmtRead())) );
00323 if( p.lgEOL() )
00324 {
00325 p.PrintLine( ioQQQ );
00326 fprintf( ioQQQ, " An observed column density MUST be entered. Sorry.\n" );
00327 cdEXIT(EXIT_FAILURE);
00328 }
00329
00330 realnum err = realnum(p.FFmtRead());
00331 if( err <= 0.0 )
00332 {
00333
00334 err = realnum(DEFERR);
00335 }
00336
00337
00338 if( p.nMatch( "<" ) )
00339 {
00340
00341 err = -err;
00342 }
00343 optimize.ColDen_error.push_back(err);
00344
00345 p.getline();
00346 if( p.m_lgEOF )
00347 {
00348 fprintf( ioQQQ, " Hit EOF while reading column density list; use END to end list.\n" );
00349 cdEXIT(EXIT_FAILURE);
00350 }
00351
00352 if( p.strcmp( "END" ) == 0 )
00353 {
00354 p.m_lgEOF = true;
00355 }
00356 }
00357
00358 if( trace.lgTrace && optimize.lgTrOpt )
00359 {
00360 fprintf( ioQQQ, "%ld columns were entered, they were;\n",
00361 (long int) optimize.ColDen_Obs.size() );
00362 for( long i=0; i < long(optimize.ColDen_Obs.size()); i++ )
00363 {
00364 fprintf( ioQQQ, " %4.4s ion=%5ld%10.2e%10.2e\n",
00365 optimize.chColDen_label[i].c_str(), optimize.ion_ColDen[i],
00366 optimize.ColDen_Obs[i], optimize.ColDen_error[i] );
00367 }
00368 }
00369 return;
00370 }
00371
00372
00373 STATIC void GetOptLineInt(Parser &p )
00374 {
00375 DEBUG_ENTRY( "GetOptLineInt()" );
00376
00377
00378 p.getline();
00379 if( p.m_lgEOF )
00380 {
00381 fprintf( ioQQQ, " Hit EOF while reading line list; use END to end list.\n" );
00382 cdEXIT(EXIT_FAILURE);
00383 }
00384
00385 while( !p.m_lgEOF )
00386 {
00387
00388 optimize.nEmergent.push_back( p.nMatch("EMER") ? 1 : 0 );
00389
00390
00391 optimize.chLineLabel.push_back( p.getCommand(4) );
00392
00393
00394 realnum wavl = realnum(p.getWaveOpt());
00395 optimize.wavelength.push_back(wavl);
00396
00397
00398 optimize.errorwave.push_back( WavlenErrorGet(wavl) );
00399
00400
00401 realnum xLineInt = realnum(p.FFmtRead());
00402 if( p.lgEOL() )
00403 {
00404 fprintf( ioQQQ, " The wavelength and relative intensity MUST be entered on this line. Sorry.\n" );
00405 fprintf( ioQQQ, " The command line is the following:\n " );
00406 p.PrintLine(ioQQQ);
00407 cdEXIT(EXIT_FAILURE);
00408 }
00409
00410 if( xLineInt <= 0. )
00411 {
00412 fprintf( ioQQQ, " An observed intensity of %.2e is not allowed. Sorry.\n",
00413 xLineInt );
00414 fprintf( ioQQQ, " The command line is the following:\n" );
00415 p.PrintLine(ioQQQ);
00416 cdEXIT(EXIT_FAILURE);
00417 }
00418 optimize.xLineInt_Obs.push_back(xLineInt);
00419
00420
00421 realnum err = realnum(p.FFmtRead());
00422
00423 if( err <= 0.0 )
00424 err = realnum(DEFERR);
00425
00426
00427 if( p.nMatch( "<" ) )
00428 {
00429
00430 err = -err;
00431 }
00432 optimize.xLineInt_error.push_back(err);
00433
00434
00435 p.getline();
00436 if( p.m_lgEOF )
00437 {
00438 fprintf( ioQQQ, " Hit EOF while reading line list for optimize command; use END to end list.\n" );
00439 cdEXIT(EXIT_FAILURE);
00440 }
00441
00442 if( p.strcmp( "END" ) == 0 )
00443 p.m_lgEOF = true;
00444 }
00445
00446 if( trace.lgTrace && trace.lgTrOptm )
00447 {
00448 fprintf( ioQQQ, "%ld lines were entered, they were:\n",
00449 (long int) optimize.xLineInt_Obs.size() );
00450
00451 for( long i=0; i < long(optimize.xLineInt_Obs.size()); i++ )
00452 {
00453 fprintf( ioQQQ, " %4.4s ", optimize.chLineLabel[i].c_str() );
00454 prt_wl( ioQQQ, optimize.wavelength[i] );
00455
00456 fprintf( ioQQQ, " %10.2e%10.2e\n",
00457 optimize.xLineInt_Obs[i],
00458 optimize.xLineInt_error[i] );
00459 }
00460 }
00461 return;
00462 }
00463
00464
00465 STATIC void GetOptTemp(Parser &p )
00466 {
00467 DEBUG_ENTRY( "GetOptTemp()" );
00468
00469
00470 p.getline();
00471 if( p.m_lgEOF )
00472 {
00473 fprintf( ioQQQ, " Hit EOF while reading line list; use END to end list.\n" );
00474 cdEXIT(EXIT_FAILURE);
00475 }
00476
00477
00478 while( !p.m_lgEOF )
00479 {
00480
00481 optimize.chTempLab.push_back( p.getCommand(4) );
00482
00483
00484 optimize.ionTemp.push_back( nint(p.FFmtRead()) );
00485
00486
00487 realnum temp_obs = realnum(p.FFmtRead());
00488 if( p.lgEOL() )
00489 {
00490 fprintf( ioQQQ, " The ion stage and temperature MUST be entered on this line. Sorry.\n" );
00491 fprintf( ioQQQ, " The command line is the following:\n " );
00492 p.PrintLine(ioQQQ);
00493 cdEXIT(EXIT_FAILURE);
00494 }
00495
00496 if( temp_obs <= 10. )
00497 temp_obs = realnum(pow( 10., temp_obs ) );
00498 optimize.temp_obs.push_back(temp_obs);
00499
00500
00501 realnum temp_error = realnum(p.FFmtRead());
00502
00503 if( temp_error <= 0.f )
00504 temp_error = realnum(DEFERR);
00505
00506 if( p.nMatch( "<" ) )
00507 temp_error = -temp_error;
00508 optimize.temp_error.push_back(temp_error);
00509
00510
00511
00512
00513
00514
00515
00516
00517 if( p.nMatch( "VOLUME" ) )
00518 optimize.chTempWeight.push_back("volume");
00519 else
00520 optimize.chTempWeight.push_back("radius");
00521
00522
00523 p.getline();
00524 if( p.m_lgEOF )
00525 {
00526 fprintf( ioQQQ, " Hit EOF while reading line list for optimize command; use END to end list.\n" );
00527 cdEXIT(EXIT_FAILURE);
00528 }
00529
00530 if( p.strcmp( "END" ) == 0 )
00531 p.m_lgEOF = true;
00532 }
00533
00534 if( trace.lgTrace && trace.lgTrOptm )
00535 {
00536 fprintf( ioQQQ, "%ld temperatures were entered, they were;\n",
00537 (long int) optimize.temp_obs.size() );
00538
00539 for( long i=0; i < long(optimize.temp_obs.size()); i++ )
00540 {
00541 fprintf( ioQQQ, " %4.4s ", optimize.chTempLab[i].c_str() );
00542 fprintf( ioQQQ, " %li " , optimize.ionTemp[i] );
00543
00544 fprintf( ioQQQ, " %.2e %.2e\n",
00545 optimize.temp_obs[i],
00546 optimize.temp_error[i] );
00547 }
00548 }
00549 return;
00550 }