00001
00002
00003
00004 #include "cddefines.h"
00005 #include "optimize.h"
00006 #include "abund.h"
00007 #include "dense.h"
00008 #include "input.h"
00009 #include "iso.h"
00010 #include "hmi.h"
00011 #include "mole.h"
00012 #include "elementnames.h"
00013 #include "parse.h"
00014
00015
00016
00017 void ParseElement(char *chCard )
00018 {
00019
00020
00021 static bool lgFirst = true;
00022 static long levels[NISO][LIMELM];
00023 char chCap[INPUT_LINE_LENGTH];
00024 bool lgEOL,
00025 lgEnd,
00026 lgHIT;
00027
00028 bool lgForceLog=false, lgForceLinear=false;
00029
00030 long int i,
00031 nelem,
00032 j,
00033 k;
00034 double param;
00035
00036 DEBUG_ENTRY( "ParseElement()" );
00037
00038
00039 if( lgFirst )
00040 {
00041 lgFirst = false;
00042 for( i=0; i<NISO; ++i )
00043 {
00044 for( nelem=0; nelem<LIMELM; ++nelem )
00045 {
00046 levels[i][nelem] = 0;
00047 }
00048 }
00049 }
00050
00051 abund.lgAbnSolar = false;
00052
00053
00054 if( nMatch("READ",chCard) )
00055 {
00056 abund.npSolar = 0;
00057 for( i=1; i < LIMELM; i++ )
00058 {
00059 input_readarray(chCard,&lgEnd);
00060
00061 if( lgEnd )
00062 {
00063 fprintf( ioQQQ, " Hit EOF while reading element list; use END to end list.\n" );
00064 puts( "[Stop in ParseElement]" );
00065 cdEXIT(EXIT_FAILURE);
00066 }
00067
00068 strcpy( chCap, chCard );
00069 caps(chCap);
00070
00071
00072 if( strncmp(chCap,"END",3) == 0 )
00073 {
00074
00075 DEBUG_EXIT( "ParseElement()" );
00076 return;
00077 }
00078
00079 j = 1;
00080 lgHIT = false;
00081 while( j < LIMELM && !lgHIT )
00082 {
00083 j += 1;
00084
00085 if( strncmp( chCap,elementnames.chElementNameShort[j-1] , 4) == 0 )
00086 {
00087 abund.npSolar += 1;
00088 abund.ipSolar[abund.npSolar-1] = j;
00089 lgHIT = true;
00090 }
00091 }
00092
00093 if( !lgHIT )
00094 {
00095 fprintf( ioQQQ, "%80.80s\n", chCard );
00096 fprintf( ioQQQ, " Sorry, but I did not recognize element name on this line.\n" );
00097 fprintf( ioQQQ, " Here is the list of names I recognize.\n" );
00098 fprintf( ioQQQ, " " );
00099
00100 for( k=2; k <= LIMELM; k++ )
00101 {
00102 fprintf( ioQQQ, "%4.4s\n", elementnames.chElementNameShort[k-1] );
00103 }
00104
00105 puts( "[Stop in ParseElement]" );
00106 cdEXIT(EXIT_FAILURE);
00107 }
00108 }
00109
00110
00111 input_readarray(chCard,&lgEnd);
00112 strcpy( chCap, chCard );
00113 caps(chCap);
00114
00115 if( strncmp(chCap,"END",3) == 0 )
00116 {
00117
00118 DEBUG_EXIT( "ParseElement()" );
00119 return;
00120 }
00121
00122 else
00123 {
00124 fprintf( ioQQQ, " Too many elements were entered.\n" );
00125 fprintf( ioQQQ, " I only know about%3d elements.\n",
00126 LIMELM );
00127 fprintf( ioQQQ, " Sorry.\n" );
00128 puts( "[Stop in ParseElement]" );
00129 cdEXIT(EXIT_FAILURE);
00130 }
00131 }
00132
00133
00134 i = 0;
00135
00136
00137 nelem = GetElem(chCard );
00138
00139 if( nelem < 0 )
00140 {
00141 fprintf( ioQQQ,
00142 " ParseElement did not find an element on the following line:\n" );
00143 fprintf( ioQQQ,
00144 "%80.80s\n", chCard );
00145 puts( "[Stop in ParseElement]" );
00146 cdEXIT(EXIT_FAILURE);
00147 }
00148
00149 ASSERT( nelem>=0 && nelem < LIMELM );
00150
00151
00152
00153 if( nMatch(" LOG",chCard) )
00154 {
00155 lgForceLog = true;
00156 }
00157 else if( nMatch("LINE",chCard) )
00158 {
00159 lgForceLinear = true;
00160 }
00161
00162 i = 4;
00163 param = FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00164
00165 if( nMatch("SCAL",chCard) )
00166 {
00167
00168 if( lgEOL )
00169 {
00170 fprintf( ioQQQ, " There must be a number on this line.\n" );
00171 fprintf( ioQQQ, " Sorry.\n" );
00172 puts( "[Stop in ParseElement]" );
00173 cdEXIT(EXIT_FAILURE);
00174 }
00175
00176 else
00177 {
00178
00179 if( (lgForceLog || param <= 0.) && !lgForceLinear )
00180 {
00181
00182 param = pow(10.,param);
00183 }
00184 abund.ScaleElement[nelem] = (float)param;
00185 }
00186 }
00187
00188 else if( nMatch("ABUN",chCard) )
00189 {
00190
00191 if( lgEOL )
00192 {
00193 fprintf( ioQQQ, " There must be a number on this line.\n" );
00194 fprintf( ioQQQ, " Sorry.\n" );
00195 puts( "[Stop in ParseElement]" );
00196 cdEXIT(EXIT_FAILURE);
00197 }
00198
00199 else
00200 {
00201 if( lgForceLinear )
00202 {
00203 abund.solar[nelem] = (float)param;
00204 }
00205 else
00206 {
00207 abund.solar[nelem] = (float)pow(10.,param);
00208 }
00209
00210 if( abund.solar[nelem]> 1. )
00211 {
00212 fprintf( ioQQQ,
00213 " Please check the abundance of this element. It seems high to me.\n" );
00214 }
00215 }
00216 }
00217
00218 else if( nMatch(" OFF",chCard) )
00219 {
00220
00221 dense.lgElmtOn[nelem] = false;
00222
00223 dense.IonHigh[nelem] = -1;
00224 dense.lgSetIoniz[nelem] = false;
00225
00226 if( nelem > ipHELIUM )
00227 {
00228 levels[ipHYDROGEN][nelem] = iso.numLevels_max[ipH_LIKE][nelem];
00229 levels[ipHELIUM][nelem] = iso.numLevels_max[ipHE_LIKE][nelem];
00230 iso.numLevels_max[ipH_LIKE][nelem] = 0;
00231 iso.numLevels_max[ipHE_LIKE][nelem] = 0;
00232 }
00233
00234 if( nelem == 0 )
00235 {
00236 fprintf( ioQQQ, " It is not possible to turn hydrogen off.\n" );
00237 fprintf( ioQQQ, " Sorry.\n" );
00238 puts( "[Stop in ParseElement]" );
00239 cdEXIT(EXIT_FAILURE);
00240 }
00241 }
00242
00243
00244 else if( nMatch("IONI",chCard) )
00245 {
00246 bool lgLogSet = false;
00247 int ion;
00248 int ihi , low;
00249
00250 if( !dense.lgElmtOn[nelem] )
00251 {
00252 fprintf(ioQQQ,"Sorry, you cannot set the ionization of %s since it has been turned off.\n" ,
00253 elementnames.chElementName[nelem] );
00254 puts( "[Stop in ParseElement]" );
00255 cdEXIT(EXIT_FAILURE);
00256 }
00257
00258 i = 4;
00259
00260 dense.lgSetIoniz[nelem] = true;
00261 ion = 0;
00262 while( ion<nelem+2 )
00263 {
00264
00265
00266
00267 dense.SetIoniz[nelem][ion] = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00268
00269 if( lgEOL )
00270 break;
00271
00272
00273 if( dense.SetIoniz[nelem][ion] < 0. )
00274 lgLogSet = true;
00275 ++ion;
00276 }
00277
00278
00279 for( i=ion; i<nelem+2; ++i )
00280 {
00281 dense.SetIoniz[nelem][i] = 0.;
00282 }
00283
00284
00285 if( lgLogSet )
00286 {
00287 for( i=0; i<ion; ++i )
00288 {
00289 dense.SetIoniz[nelem][i] = (float)pow(10.f, dense.SetIoniz[nelem][i]);
00290 }
00291 }
00292
00293
00294 low = 0;
00295 while( dense.SetIoniz[nelem][low]==0 && low < nelem+1 )
00296 ++low;
00297 ihi = nelem+1;
00298 while( dense.SetIoniz[nelem][ihi]==0 && ihi > low )
00299 --ihi;
00300
00301 if( ihi==low && dense.SetIoniz[nelem][ihi]==0 )
00302 {
00303 fprintf(ioQQQ," element ionization command has all zero ionization fractions. This is not possible.\n Sorry\n");
00304 puts( "[Stop in ParseElement]" );
00305 cdEXIT(EXIT_FAILURE);
00306 }
00307 for(i=low; i<=ihi; ++i )
00308 {
00309 if( dense.SetIoniz[nelem][i]==0 )
00310 {
00311 fprintf(ioQQQ," element abundance command has zero abundance between positive values. This is not possible.\n Sorry\n");
00312 puts( "[Stop in ParseElement]" );
00313 cdEXIT(EXIT_FAILURE);
00314 }
00315 }
00316
00317
00318
00319
00320
00321
00322 if( nelem==ipHYDROGEN && !hmi.lgNoH2Mole )
00323 {
00324 fprintf(ioQQQ," Hydrogen ionization has been set, H chemistry is disabled.\n");
00325
00326 hmi.lgNoH2Mole = true;
00327 }
00328
00329
00330
00331 if( mole.lgElem_in_chemistry[nelem] && !co.lgNoCOMole )
00332 {
00333 fprintf(ioQQQ," The ionization of an element in the CO chemistry network has been set, CO chemistry is disabled.\n");
00334
00335 co.lgNoCOMole = true;
00336 }
00337 }
00338
00339
00340 else if( nMatch(" ON ",chCard) )
00341 {
00342
00343 dense.lgElmtOn[nelem] = true;
00344
00345 if( levels[ipHYDROGEN][nelem] )
00346 {
00347 iso.numLevels_max[ipH_LIKE][nelem] = levels[ipHYDROGEN][nelem];
00348 iso.numLevels_max[ipHE_LIKE][nelem] = levels[ipHELIUM][nelem];
00349 }
00350 }
00351
00352 else if( nMatch("TABL",chCard) )
00353 {
00354
00355
00356 abund.lgAbunTabl[nelem] = true;
00357
00358
00359 abund.lgAbTaON = true;
00360 if( nMatch("DEPT",chCard) )
00361 {
00362 abund.lgAbTaDepth[nelem] = true;
00363 }
00364 else
00365 {
00366 abund.lgAbTaDepth[nelem] = false;
00367 }
00368
00369
00370 if( nelem == 0 )
00371 {
00372 fprintf( ioQQQ, " cannot change abundance of hydrogen.\n" );
00373 fprintf( ioQQQ, " Sorry.\n" );
00374 puts( "[Stop in ParseElement]" );
00375 cdEXIT(EXIT_FAILURE);
00376 }
00377
00378 input_readarray(chCard,&lgEnd);
00379 i = 1;
00380 abund.AbTabRad[0][nelem] = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00381 abund.AbTabFac[0][nelem] = (float)FFmtRead(chCard,&i,INPUT_LINE_LENGTH,&lgEOL);
00382
00383 if( lgEOL )
00384 {
00385 fprintf( ioQQQ, " no pairs entered - cant interpolate\n" );
00386 puts( "[Stop in ParseElement]" );
00387 cdEXIT(EXIT_FAILURE);
00388 }
00389
00390 abund.nAbunTabl = 2;
00391 lgEnd = false;
00392
00393 while( !lgEnd && abund.nAbunTabl < LIMTABD )
00394 {
00395 input_readarray(chCard,&lgEnd);
00396 if( !lgEnd )
00397 {
00398
00399 cap4(chCap , chCard);
00400 if( strncmp(chCap,"END",3) == 0 )
00401 lgEnd = true;
00402 }
00403
00404
00405 if( !lgEnd )
00406 {
00407 i = 1;
00408 abund.AbTabRad[abund.nAbunTabl-1][nelem] =
00409 (float)FFmtRead(chCard ,&i,INPUT_LINE_LENGTH,&lgEOL);
00410
00411 abund.AbTabFac[abund.nAbunTabl-1][nelem] =
00412 (float)FFmtRead(chCard ,&i,INPUT_LINE_LENGTH,&lgEOL);
00413 abund.nAbunTabl += 1;
00414 }
00415 }
00416
00417 abund.nAbunTabl -= 1;
00418
00419
00420 for( i=1; i < abund.nAbunTabl; i++ )
00421 {
00422
00423 if( abund.AbTabRad[i][nelem] <= abund.AbTabRad[i-1][nelem] )
00424 {
00425 fprintf( ioQQQ, " abun radii must be in increasing order\n" );
00426 puts( "[Stop in ParseElement]" );
00427 cdEXIT(EXIT_FAILURE);
00428 }
00429 }
00430 }
00431
00432 else
00433 {
00434 fprintf( ioQQQ, " There must be a keyword on this line.\n" );
00435 fprintf( ioQQQ, " The keys I know about are TABLE, SCALE, _OFF, _ON_, IONIZATION, and ABUNDANCE.\n" );
00436 fprintf( ioQQQ, " Sorry.\n" );
00437 puts( "[Stop in ParseElement]" );
00438 cdEXIT(EXIT_FAILURE);
00439 }
00440
00441
00442 if( optimize.lgVarOn )
00443 {
00444 optimize.nvarxt[optimize.nparm] = 1;
00445
00446 optimize.nvfpnt[optimize.nparm] = input.nRead;
00447
00448 if( nMatch("SCAL",chCard) )
00449 {
00450
00451 sprintf( optimize.chVarFmt[optimize.nparm], "ELEMENT %4.4s SCALE %%f LOG",
00452 elementnames.chElementNameShort[nelem] );
00453
00454
00455 optimize.vparm[0][optimize.nparm] = (float)log10(param);
00456 optimize.vincr[optimize.nparm] = 0.2f;
00457 }
00458
00459 else if( nMatch("ABUN",chCard) )
00460 {
00461
00462 sprintf( optimize.chVarFmt[optimize.nparm], "ELEMENT %4.4s ABUND %%f LOG ",
00463 elementnames.chElementNameShort[nelem] );
00464
00465
00466 optimize.vparm[0][optimize.nparm] = (float)param;
00467 optimize.vincr[optimize.nparm] = 0.2f;
00468 }
00469 ++optimize.nparm;
00470 }
00471
00472 DEBUG_EXIT( "ParseElement()" );
00473 return;
00474 }
00475
00476
00477