57 if( strcmp(
"H2 ",
label.c_str() ) == 0 )
98 H2_SaveLine[(*Hi).n()][(*Hi).v()][(*Hi).J()][(*Lo).n()][(*Lo).v()][(*Lo).J()] = 0.;
102 H2_SaveLine[(*Hi).n()][(*Hi).v()][(*Hi).J()][(*Lo).n()][(*Lo).v()][(*Lo).J()] +=
116 save.whichDiatomToPrint[
save.nsave] = &(*this);
122 strcpy(
save.chSave[
save.nsave],
"H2cl" );
129 "H2 vibration state",0.0);
133 "H2 rotation state",0.0);
140 sprintf( chHeader,
"#vib\trot\tcolumn density\n" );
146 sprintf( chHeader,
"#vib\trot\tEner(K)\tcolden\tcolden/stat wght\tLTE colden\tLTE colden/stat wght\n" );
149 else if( p.
nMatch(
"COOL") )
152 strcpy(
save.chSave[
save.nsave],
"H2co" );
154 "#H2 depth\ttot cool\tTH Sol\tBig Sol\tTH pht dis\tpht dis\tTH Xcool\tXcool \n" );
157 else if( p.
nMatch(
"CREA") )
160 fprintf(
ioQQQ,
" This command has been superseded by the \"creation\" option of the \"save chemistry rates\" command.\n" );
161 fprintf(
ioQQQ,
" Sorry.\n" );
164 else if( p.
nMatch(
"DEST") )
167 fprintf(
ioQQQ,
" This command has been superseded by the \"destruction\" option of the \"save chemistry rates\" command.\n" );
168 fprintf(
ioQQQ,
" Sorry.\n" );
172 else if( p.
nMatch(
"HEAT") )
175 strcpy(
save.chSave[
save.nsave],
"H2he" );
177 "#H2 depth\ttot Heat\tHeat(big)\tHeat(TH85)\tDissoc(Big)\tDissoc(TH85) \n" );
180 else if( p.
nMatch(
"LEVE") )
183 strcpy(
save.chSave[
save.nsave],
"H2le" );
185 "#H2 v\tJ\tenergy(wn)\tstat wght\tSum As" );
191 strcat( chHeader , chHoldit );
193 strcat( chHeader ,
"\n" );
196 else if( p.
nMatch(
"LINE") )
199 strcpy(
save.chSave[
save.nsave],
"H2ln" );
201 "#H2 line\tEhi\tVhi\tJhi\tElo\tVlo\tJlo\twl(mic)\twl(lab)\tlog L or I\tI/Inorm\tExcit(hi, K)\tg_u h nu * Aul\n" );
208 "faintest line to save",1e-4);
223 else if( p.
nMatch(
"GROU") )
231 "electronic levels for output",1.0);
236 else if( p.
nMatch(
" PDR") )
239 strcpy(
save.chSave[
save.nsave],
"H2pd" );
240 sprintf( chHeader,
"#H2 creation, destruction. \n" );
242 else if( p.
nMatch(
"POPU") )
245 strcpy(
save.chSave[
save.nsave],
"H2po" );
252 "highest H2 save vibration state",0.0);
256 "highest H2 save rotation state",0.0);
262 sprintf( chHeader,
"#depth\torth\tpar\te=1 rel pop\te=2 rel pop\tv,J rel pops\n" );
273 sprintf( chHeader,
"#vib\trot\tpops\n" );
279 sprintf( chHeader,
"#vib\trot\ts\tenergy(wn)\tpops/H2\told/H2\tpops/g/H2\tdep coef\tFin(Col)\tFout(col)\tRCout\tRRout\tRCin\tRRin\n" );
284 else if( p.
nMatch(
"RATE") )
287 strcpy(
save.chSave[
save.nsave],
"H2ra" );
289 "#depth\tN(H2)\tN(H2)/u(H2)\tA_V(star)\tn(Eval)"
290 "\tH2/Htot\trenorm\tfrm grn\tfrmH-\tdstTH85\tBD96\tELWERT\tBigH2\telec->H2g\telec->H2s"
291 "\tG(TH85)\tG(DB96)\tCR\tEleclife\tShield(BD96)\tShield(H2)\tBigh2/G0(spc)\ttot dest"
292 "\tHeatH2Dish_TH85\tHeatH2Dexc_TH85\tHeatDish_BigH2\tHeatDexc_BigH2\thtot\n" );
294 else if( p.
nMatch(
"SOLO") )
297 strcpy(
save.chSave[
save.nsave],
"H2so" );
299 "#depth\tSol tot\tpump/dissoc\tpump/dissoc BigH2\tavH2g\tavH2s\tH2g chem/big H2\tH2s chem/big H2\tfrac H2g BigH2\tfrac H2s BigH2\teHi\tvHi\tJHi\tvLo\tJLo\tfrac\twl(A)\n" );
301 else if( p.
nMatch(
"SPEC") )
304 strcpy(
save.chSave[
save.nsave],
"H2sp" );
306 "#depth\tspecial\n" );
308 else if( p.
nMatch(
"TEMP") )
311 strcpy(
save.chSave[
save.nsave],
"H2te" );
313 "#depth\tH2/H\tn(1/0)\tn(ortho/para)\tT(1/0)\tT(2/0)\tT(3/0)\tT(3/1)\tT(4/0)\tT(kin)\tT(21cm)\tT_sum(1/0)\tT_sum(2/0)\tT_sum(3/0)\tT_sum(3/1)\tT_sum(4/0) \n");
315 else if( p.
nMatch(
"THER") )
318 strcpy(
save.chSave[
save.nsave],
"H2th" );
320 "#depth\tH2/H\tn(1/0)\tn(ortho/para)\tT(1/0)\tT(2/0)\tT(3/0)\tT(3/1)\tT(4/0)\tT(kin)\tT(21cm)\tT_sum(1/0)\tT_sum(2/0)\tT_sum(3/0)\tT_sum(3/1)\tT_sum(4/0) \n");
325 " There must be a second key; they are RATE, LINE, COOL, COLUMN, _PDR, SOLOmon, TEMP, and POPUlations\n" );
341 fprintf(
ioQQQ,
" %s density ",
label.c_str() );
344 fprintf(
ioQQQ,
" orth/par");
347 fprintf(
ioQQQ,
" v0 J=0,3");
353 fprintf(
ioQQQ,
" TOTv=0,3");
358 fprintf(
ioQQQ,
"\n");
371 fprintf(
ioQQQ,
" %s departure coefficients\n",
label.c_str() );
374 fprintf(
ioQQQ,
"%li electronic\n", iElec );
375 for(
long iVib=0; iVib<=
nVib_hi[iElec]; ++iVib )
377 for(
long iRot=0; iRot<
Jlowest[iElec]; ++iRot )
378 fprintf(
ioQQQ,
" -----" );
379 for(
long iRot=
Jlowest[iElec]; iRot<=
nRot_hi[iElec][iVib]; ++iRot )
384 fprintf(
ioQQQ,
"\n" );
386 fprintf(
ioQQQ,
"\n" );
409 fprintf( ioMEAN,
" H2 total ");
412 fprintf( ioMEAN,
" H2 ortho ");
415 fprintf( ioMEAN,
" para");
419 fprintf( ioMEAN,
" v0 J=0,3");
432 const char* cdDATAFILE[
N_ELEC] =
436 "transprob_C_plus.dat",
437 "transprob_C_minus.dat",
438 "transprob_B_primed.dat",
439 "transprob_D_plus.dat",
440 "transprob_D_minus.dat"
444 long int i, n1, n2, n3;
445 long int iVibHi , iVibLo , iRotHi , iRotLo , iElecHi , iElecLo;
452 strcpy( chPath,
path.c_str() );
453 strcat( chPath,
input.chDelimiter );
454 strcat( chPath, cdDATAFILE[nelec] );
460 fprintf(
ioQQQ,
" H2_ReadTransprob could not read first line of %s\n", cdDATAFILE[nelec]);
465 n1 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
466 n2 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
467 n3 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
471 if( ( n1 != 2 ) || ( n2 != 4 ) || ( n3 != 29 ) )
474 " H2_ReadTransprob: the version of %s is not the current version.\n", cdDATAFILE[nelec] );
476 " I expected to find the number 2 4 29 and got %li %li %li instead.\n" ,
478 fprintf(
ioQQQ,
"Here is the line image:\n==%s==\n", chLine );
483 while(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) != NULL )
488 if( chLine[0]==
'\n' || chLine[0]==
'\0' || chLine[0]==
' ' )
492 int n = sscanf(chLine,
"%li\t%li\t%li\t%li\t%li\t%li\t%le",
493 &iElecHi , &iVibHi ,&iRotHi , &iElecLo , &iVibLo , &iRotLo , &Aul );
495 ASSERT( iElecHi == nelec );
500 if( iVibHi <=
nVib_hi[iElecHi] &&
502 iRotHi <=
nRot_hi[iElecHi][iVibHi] &&
503 iRotLo <=
nRot_hi[iElecLo][iVibLo])
507 double ener =
states[ipHi].energy().WN() -
states[ipLo].energy().WN();
511 trns[lineIndex].AddLine2Stack();
521 fprintf(
ioQQQ,
"negative energy H2 transition\t%li\t%li\t%li\t%li\t%.2e\t%.2e\n",
522 iVibHi,iVibLo,iRotHi,iRotLo,Aul,ener);
535void H2_Read_Cosmicray_distribution(
void)
542 long int i, n1, n2, n3, iVib , iRot;
550 strcpy( chPath, path.c_str() );
551 strcat( chPath,
input.chDelimiter );
552 strcat( chPath,
"H2_CosmicRay_collision.dat" );
558 fprintf(
ioQQQ,
" H2_Read_Cosmicray_distribution could not read first line of %s\n",
"H2_Cosmic_collision.dat");
564 n1 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
565 n2 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
566 n3 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
570 if( ( n1 != 1 ) || ( n2 != 21 ) || ( n3 != 3 ) )
573 " H2_Read_Cosmicray_distribution: the version of %s is not the current version.\n",
"H2_Cosmic_collision.dat" );
575 " I expected to find the number 1 21 3 and got %li %li %li instead.\n" ,
577 fprintf(
ioQQQ,
"Here is the line image:\n==%s==\n", chLine );
585 while( chLine[0]==
'#' )
599 sscanf(chLine,
"%li\t%li\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf",
600 &iVib ,&j_minus_ji , &a[0],&a[1],&a[2],&a[3],&a[4],&a[5],&a[6],&a[7],&a[8],&a[9]
609 ASSERT( j_minus_ji == -2 || j_minus_ji == +2 || j_minus_ji == 0 );
610 ASSERT( neut_frac < CR_X );
613 j_minus_ji = 1 + j_minus_ji/2;
614 ASSERT( j_minus_ji>=0 && j_minus_ji<=2 );
617 for( iRot=0; iRot<CR_J; ++iRot )
619 cr_rate[neut_frac][iVib][iRot][j_minus_ji] = (
realnum)a[iRot];
621 if( lgH2_NOISECOSMIC )
626 for( iRot=0; iRot<CR_J; ++iRot )
628 cr_rate[neut_frac][iVib][iRot][j_minus_ji] *= (
realnum)pow(10.,(
double)r);
634 fprintf(
ioQQQ,
"cr rate\t%li\t%li", iVib , j_minus_ji );
635 for( iRot=0; iRot<CR_J; ++iRot )
637 fprintf(
ioQQQ,
"\t%.3e", cr_rate[neut_frac][iVib][iRot][j_minus_ji] );
639 fprintf(
ioQQQ,
"\n" );
645 while( chLine[0]==
'#' )
685 vector<level_tmp> levels;
686 levels.resize( n.size() );
687 ASSERT( levels.size() > 0 );
688 for(
unsigned i = 0; i < n.size(); ++i )
693 levels[i].eWN = eWN[i];
697 sort( levels.begin(), levels.end() );
700 for( vector<level_tmp>::iterator lev = levels.begin(); lev != levels.end(); ++lev )
703 long i =
states.size() - 1;
707 states[i].energy().set( lev->eWN,
"cm^-1" );
731 nRot_hi[ (*st).n() ][ (*st).v() ] =
MAX2(
nRot_hi[ (*st).n() ][ (*st).v() ], (*st).J() );
740 const char* cdDATAFILE[
N_ELEC] =
745 "energy_C_minus.dat",
746 "energy_B_primed.dat",
752 strcpy( chPath,
path.c_str() );
753 strcat( chPath,
input.chDelimiter );
754 strcat( chPath, cdDATAFILE[nelec] );
758 long int i, n1, n2, n3;
764 fprintf(
ioQQQ,
" H2_ReadEnergies could not read first line of %s\n", cdDATAFILE[nelec]);
769 n1 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
770 n2 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
771 n3 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
775 if( ( n1 != 2 ) || ( n2 != 4 ) || ( n3 != 29 ) )
778 " H2_ReadEnergies: the version of %s is not the current version.\n", cdDATAFILE[nelec] );
780 " I expected to find the number 2 4 29 and got %li %li %li instead.\n" ,
782 fprintf(
ioQQQ,
"Here is the line image:\n==%s==\n", chLine );
791 while(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) != NULL )
796 if( chLine[0]==
'\n' || chLine[0]==
'\0' || chLine[0]==
' ' )
800 int nReads = sscanf(chLine,
"%li\t%li\t%le", &iVib, &iRot, &energyWN );
804 ASSERT( energyWN > 0. || (nelec==0 && iVib==0 && iRot==0 ) );
806 n.push_back( nelec );
809 eWN.push_back( energyWN );
831 const char* cdDATAFILE =
"energy_dissoc.dat";
834 long int i, n1, n2, n3;
841 strcpy( chPath,
path.c_str() );
842 strcat( chPath,
input.chDelimiter );
843 strcat( chPath, cdDATAFILE );
849 fprintf(
ioQQQ,
" H2_ReadDissocEnergies could not read first line of %s\n", cdDATAFILE );
854 n1 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
855 n2 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
856 n3 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
860 if( ( n1 != 2 ) || ( n2 != 4 ) || ( n3 != 29 ) )
863 " H2_ReadDissocEnergies: the version of %s is not the current version.\n", cdDATAFILE );
865 " I expected to find the number 2 4 29 and got %li %li %li instead.\n" ,
867 fprintf(
ioQQQ,
"Here is the line image:\n==%s==\n", chLine );
871 while(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) != NULL )
876 if( chLine[0]==
'\n' || chLine[0]==
'\0' || chLine[0]==
' ' )
880 int n = sscanf(chLine,
"%li\t%le", &iElec, &energyWN );
895 const char* cdDATAFILE[
N_ELEC] =
899 "dissprob_C_plus.dat",
900 "dissprob_C_minus.dat",
901 "dissprob_B_primed.dat",
902 "dissprob_D_plus.dat",
903 "dissprob_D_minus.dat"
914 strcpy( chPath,
path.c_str() );
915 strcat( chPath,
input.chDelimiter );
916 strcat( chPath, cdDATAFILE[nelec] );
922 fprintf(
ioQQQ,
" H2_ReadDissprob could not read first line of %s\n", cdDATAFILE[nelec]);
927 long n1 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
928 long n2 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
929 long n3 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
933 if( ( n1 != 3 ) || ( n2 != 2 ) || ( n3 != 11 ) )
936 " H2_ReadDissprob: the version of %s is not the current version.\n", cdDATAFILE[nelec] );
938 " I expected to find the number 3 2 11 and got %li %li %li instead.\n" ,
940 fprintf(
ioQQQ,
"Here is the line image:\n==%s==\n", chLine );
944 while(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) != NULL )
949 if( chLine[0]==
'\n' || chLine[0]==
'\0' || chLine[0]==
' ' )
955 sscanf(chLine,
"%li\t%li\t%le\t%le",
970 ( iRot >
nRot_hi[nelec][iVib] ) )
988 long int i, n1, n2, n3, iVib , iRot;
992 const bool lgH2HMINUS_PRT =
false;
998 strcpy( chPath,
path.c_str() );
999 strcat( chPath,
input.chDelimiter );
1000 strcat( chPath,
"hminus_deposit.dat" );
1004 if(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) == NULL )
1006 fprintf(
ioQQQ,
" H2_Read_hminus_distribution could not read first line of %s\n", chPath );
1012 n1 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
1013 n2 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
1014 n3 = (long)
FFmtRead(chLine,&i,
sizeof(chLine),&lgEOL);
1018 if( ( n1 != 2 ) || ( n2 != 10 ) || ( n3 != 17 ) )
1021 " H2_Read_hminus_distribution: the version of %s is not the current version.\n", chPath );
1023 " I expected to find the number 2 10 17 and got %li %li %li instead.\n" ,
1025 fprintf(
ioQQQ,
"Here is the line image:\n==%s==\n", chLine );
1030 if(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) == NULL )
1033 while( chLine[0]==
'#' )
1035 if(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) == NULL )
1046 sscanf(chLine,
"%li\t%li\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf\t%lf",
1047 &iVib ,&iRot , &ener, &a[0],&a[1],&a[2] , &a[3],&a[4],&a[5] ,&a[6]
1057 if( lgH2HMINUS_PRT )
1058 fprintf(
ioQQQ,
"hminusss\t%li\t%li", iVib , iRot );
1063 if( lgH2HMINUS_PRT )
1066 if( lgH2HMINUS_PRT )
1067 fprintf(
ioQQQ,
"\n" );
1069 if(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) == NULL )
1071 while( chLine[0]==
'#' )
1073 if(
read_whole_line( chLine , (
int)
sizeof(chLine) , ioDATA ) == NULL )
1079 if( lgH2HMINUS_PRT )
1082 fprintf(
ioQQQ,
" total H- formation rate ");
1086 fprintf(
ioQQQ,
"\t%.3e" , sumrate[i]);
1088 fprintf(
ioQQQ,
"\n" );
1092 for( iVib=0; iVib<=
nVib_hi[0]; ++iVib )
1103 if( lgH2HMINUS_PRT )
1106 fprintf(
ioQQQ,
" H- distribution function ");
1107 for( iVib=0; iVib<=
nVib_hi[0]; ++iVib )
1111 fprintf(
ioQQQ,
"%li\t%li", iVib , iRot );
1116 fprintf(
ioQQQ,
"\n" );
1139 " H2_Punch_line_data ALL option not implemented in H2_Punch_line_data yet 1\n" );
1144 bool lgPrint =
false;
1145 fprintf( ioPUN,
"#Eu\tVu\tJu\tEl\tVl\tJl\tWL\tgl\tgu\tgf\tA\tCS\tn(crt)\n" );
1149 if( (*tr).ipCont() <= 0 )
1151 (*tr).Coll().col_str() = 0.;
1155 fprintf(ioPUN,
"%2li\t%2li\t%2li\t%2li\t%2li\t%2li\t",
1156 (*Hi).n(), (*Hi).v(), (*Hi).J(),
1157 (*Lo).n(), (*Lo).v(), (*Lo).J() );
1161 fprintf( ioPUN ,
"\n");
1177 if( (*tr).ipCont() <= 0 )
1197 if( (*tr).ipCont() <= 0 )
1210 char chBranch[5] = {
'O',
'P',
'Q',
'R',
'S'};
1212 int ip = 2 + (iRotHi - iRotLo);
1215 fprintf(
ioQQQ,
" chMolBranch called with insane iRotHi=%li iRotLo=%li ip=%i\n",
1216 iRotHi , iRotLo , ip );
1220 return( chBranch[ip] );
1233 if( (strcmp( chJOB ,
"H2po" ) == 0) && (strcmp(chTime,
"LAST") == 0) &&
1234 (
save.punarg[ipPun][2] != 0) )
1242 long LimVib, LimRot;
1247 if(
save.punarg[ipPun][0] > 0 )
1249 LimVib = (long)
save.punarg[ipPun][0];
1257 fprintf(io,
"%i\t%i\t%.3e\tortho\n",
1261 fprintf(io,
"%i\t%i\t%.3e\tpara\n",
1265 fprintf(io,
"%i\t%i\t%.3e\ttotal\n",
1271 for( iVibHi=0; iVibHi<=LimVib; ++iVibHi )
1274 if(
save.punarg[ipPun][1] > 0 )
1276 LimRot = (long)
MIN2(
1281 LimRot =
nRot_hi[iElecHi][iVibHi];
1283 if(
save.punarg[ipPun][2] > 0 )
1289 fprintf(io,
"vib\\rot");
1291 for( i=0; i<=LimRot; ++i )
1293 fprintf(io,
"\t%li",i);
1297 fprintf(io,
"%li",iVibHi );
1298 for( iRotHi=
Jlowest[iElecHi]; iRotHi<=LimRot; ++iRotHi )
1300 fprintf(io,
"\t%.3e",
1305 else if(
save.punarg[ipPun][2] < 0 )
1308 for( iRotHi=
Jlowest[iElecHi]; iRotHi<=LimRot; ++iRotHi )
1313 const char chlgPara[2]={
'P',
'O'};
1314 const long ipHi =
ipEnergySort[iElecHi][iVibHi][iRotHi];
1317 fprintf(io,
"%li\t%li\t%c\t%.1f\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\n",
1321 chlgPara[
H2_lgOrtho[iElecHi][iVibHi][iRotHi]],
1323 states[ipHi].energy().WN(),
1353 else if( (strcmp( chJOB ,
"H2po" ) == 0) && (strcmp(chTime,
"LAST") != 0) &&
1354 (
save.punarg[ipPun][2] == 0) )
1359 fprintf(io,
"%.5e\t%.3e\t%.3e",
radius.depth_mid_zone ,
1362 fprintf(io,
"\t%.3e\t%.3e",
1366 long LimVib, LimRot;
1368 if(
save.punarg[ipPun][0] > 0 )
1370 LimVib = (long)
save.punarg[ipPun][1];
1374 LimVib =
nRot_hi[iElecHi][iVibHi];
1378 if(
save.punarg[ipPun][1] > 0 )
1380 LimRot = (long)
save.punarg[ipPun][1];
1384 LimRot =
nRot_hi[iElecHi][iVibHi];
1386 for( iVibHi = 0; iVibHi<=LimVib; ++iVibHi )
1388 fprintf(io,
"\tv=%li",iVibHi);
1389 long int LimRotVib =
MIN2( LimRot ,
nRot_hi[iElecHi][iVibHi] );
1390 for(
long iRotHi=
Jlowest[iElecHi]; iRotHi<=LimRotVib; ++iRotHi )
1392 fprintf(io,
"\t%.3e",
1401 else if( (strcmp( chJOB ,
"H2cl" ) == 0) && (strcmp(chTime,
"LAST") == 0) )
1406 long LimVib, LimRot;
1411 if(
save.punarg[ipPun][0] > 0 )
1413 LimVib = (long)
save.punarg[ipPun][0];
1421 fprintf(io,
"%i\t%i\t%.3e\tortho\n",
1425 fprintf(io,
"%i\t%i\t%.3e\tpara\n",
1430 fprintf(io,
"%i\t%i\t%.3e\ttotal\n",
1436 for( iVibHi=0; iVibHi<=LimVib; ++iVibHi )
1441 if(
save.punarg[ipPun][1] > 0 )
1443 LimRot = (long)
save.punarg[ipPun][1];
1447 LimRot =
nRot_hi[iElecHi][iVibHi];
1449 if(
save.punarg[ipPun][2] > 0 )
1455 fprintf(io,
"vib\\rot");
1457 for( i=0; i<=LimRot; ++i )
1459 fprintf(io,
"\t%li",i);
1463 fprintf(io,
"%li",iVibHi );
1464 for( iRotHi=
Jlowest[iElecHi]; iRotHi<=LimRot; ++iRotHi )
1466 fprintf(io,
"\t%.3e",
1474 for( iRotHi=
Jlowest[iElecHi]; iRotHi<=LimRot; ++iRotHi )
1476 fprintf(io,
"%li\t%li\t%.1f\t%.3e\t%.3e\t%.3e\t%.3e\n",
1492 else if( (strcmp(chJOB ,
"H2pd" ) == 0) && (strcmp(chTime,
"LAST") != 0) )
1496 fprintf(io,
"%.5e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
1503 hmi.H2_Solomon_dissoc_rate_TH85_H2g ,
1505 hmi.H2_Solomon_dissoc_rate_BD96_H2g,
1509 else if( (strcmp(chJOB ,
"H2co" ) == 0) && (strcmp(chTime,
"LAST") != 0) )
1512 fprintf(io,
"%.5e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
1518 hmi.H2_Solomon_dissoc_rate_TH85_H2g,
1523 hmi.HeatH2Dish_TH85,
1527 hmi.HeatH2Dexc_TH85,
1533 else if( (strcmp(chJOB ,
"H2le" ) == 0) && (strcmp(chTime,
"LAST") == 0) )
1545 for(
long int ipLo=0; ipLo<ipHi; ++ipLo )
1553 if( ( abs(iRotHi-iRotLo) == 2 || (iRotHi-iRotLo) == 0 ) && iVibLo <= iVibHi &&
1561 Csum[nColl] += H2cr[nColl];
1565 fprintf(io,
"%li\t%li\t%.2f\t%li\t%.3e",
1567 states[ipHi].energy().WN(),
1572 fprintf(io,
"\t%.3e",Csum[nColl]);
1577 else if( (strcmp(chJOB ,
"H2ra" ) == 0) && (strcmp(chTime,
"LAST") != 0) )
1580 double sumpop = 0. , sumlife = 0.;
1589 if( (*Lo).n() > 0 || (*Lo).v() > 0 )
1591 sumlife += (*tr).Emis().pump() * (*(*tr).Lo()).Pop();
1592 sumpop += (*(*tr).Lo()).Pop();
1600 "%.5e\t%.3e\t%.3e\t%.3e\t%li",
1609 rfield.extin_mag_V_point,
1613 "\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e",
1619 gv.rate_h2_form_grains_used_total ,
1623 hmi.H2_Solomon_dissoc_rate_TH85_H2g +
hmi.H2_Solomon_dissoc_rate_TH85_H2s,
1625 hmi.H2_Solomon_dissoc_rate_BD96_H2g +
hmi.H2_Solomon_dissoc_rate_BD96_H2s,
1627 hmi.H2_Solomon_dissoc_rate_ELWERT_H2g +
hmi.H2_Solomon_dissoc_rate_ELWERT_H2g,
1636 "\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e",
1638 hmi.UV_Cont_rel2_Habing_TH85_depth,
1640 hmi.UV_Cont_rel2_Draine_DB96_depth,
1643 sumlife/
SDIV( sumpop ) ,
1644 hmi.H2_Solomon_dissoc_rate_BD96_H2g/
SDIV(
hmi.UV_Cont_rel2_Habing_TH85_depth) ,
1647 hmi.H2_rate_destroy);
1649 "\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
1650 hmi.HeatH2Dish_TH85,
1651 hmi.HeatH2Dexc_TH85,
1657 else if( (strcmp(chJOB ,
"H2so" ) == 0) && (strcmp(chTime,
"LAST") != 0) )
1660 const int nSOL = 100;
1662 long int jlosave[nSOL] , ivlosave[nSOL],
1663 iehisave[nSOL] ,jhisave[nSOL] , ivhisave[nSOL],
1667 realnum fsave[nSOL], wlsave[nSOL];
1669 fprintf(io,
"%.5e\t%.3e",
1684 sum += (*(*tr).Lo()).Pop() * (*tr).Emis().pump();
1691 const double frac = 0.01;
1698 one = (*(*tr).Lo()).Pop() * (*tr).Emis().pump();
1699 if( one/sum > frac && nsave<nSOL)
1702 fsave[nsave] = (
realnum)(one/sum);
1703 jlosave[nsave] = (*Lo).J();
1704 ivlosave[nsave] = (*Lo).v();
1705 jhisave[nsave] = (*Hi).J();
1706 ivhisave[nsave] = (*Hi).v();
1707 iehisave[nsave] = (*Hi).n();
1708 wlsave[nsave] = (*tr).WLAng();
1731 fprintf(io,
"\t%.3f\t%.3f\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e",
1742 for(
long i=0; i<nsave; ++i )
1744 long ip = ipOrdered[i];
1746 fprintf(io,
"\t%li\t%li\t%li\t%li\t%li\t%.3f\t%.3f",
1747 iehisave[ip],ivhisave[ip],jhisave[ip],ivlosave[ip] , jlosave[ip] , fsave[ip] , wlsave[ip] );
1754 else if( (strcmp(chJOB ,
"H2te" ) == 0) && (strcmp(chTime,
"LAST") != 0) )
1757 double pop_ratio10,pop_ratio20,pop_ratio30,pop_ratio31,pop_ratio40;
1758 double T10,T20,T30,T31,T40;
1760 double T10_sum,T20_sum,T30_sum,T31_sum,T40_sum;
1761 double pop_ratio10_sum,pop_ratio20_sum,pop_ratio30_sum,pop_ratio31_sum,pop_ratio40_sum;
1764 double pop0 =
states[0].Pop();
1765 double pop1 =
states[1].Pop();
1766 double pop2 =
states[2].Pop();
1767 double pop3 =
states[3].Pop();
1768 double pop4 =
states[4].Pop();
1770 double energyK =
states[1].energy().K() -
states[0].energy().K();
1772 pop_ratio10 = pop1/
SDIV(pop0);
1778 energyK =
states[2].energy().K() -
states[0].energy().K();
1779 pop_ratio20 = pop2/
SDIV(pop0);
1785 energyK =
states[3].energy().K() -
states[0].energy().K();
1786 pop_ratio30 = pop3/
SDIV(pop0);
1792 energyK =
states[3].energy().K() -
states[1].energy().K();
1793 pop_ratio31 = pop3/
SDIV(pop1);
1799 energyK =
states[4].energy().K() -
states[0].energy().K();
1800 pop_ratio40 = pop4/
SDIV(pop0);
1809 pop_ratio10_sum = 0.;
1824 "%.5e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\n" ,
1833 T10,T20,T30,T31,T40,
1835 hyperfine.Tspin21cm,T10_sum,T20_sum,T30_sum,T31_sum,T40_sum );
1837 else if( (strcmp(chJOB ,
"H2ln" ) == 0) && (strcmp(chTime,
"LAST") == 0) )
1870 long iElecHi = (*Hi).n();
1871 long iVibHi = (*Hi).v();
1872 long iRotHi = (*Hi).J();
1873 long iElecLo = (*Lo).n();
1874 long iVibLo = (*Lo).v();
1875 long iRotLo = (*Lo).J();
1878 if(
H2_SaveLine[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo] > thresh )
1882 double wl = (*tr).WLAng()/1e4;
1883 fprintf(io,
"%li-%li %c(%li)", iVibHi, iVibLo,
chMolBranch( iRotHi, iRotLo ), iRotLo );
1884 fprintf( io,
"\t%ld\t%ld\t%ld\t%ld\t%ld\t%ld", iElecHi , iVibHi , iRotHi , iElecLo , iVibLo , iRotLo);
1886 fprintf( io,
"\t%.7f\t", wl );
1888 prt_wl( io , (*tr).WLAng() );
1890 fprintf( io,
"\t%.3f\t%.3e",
1891 log10(
MAX2(1e-37,
H2_SaveLine[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo])) +
radius.Conv2PrtInten,
1892 H2_SaveLine[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo]*renorm );
1894 fprintf( io,
"\t%.3f", (*Hi).energy().K() );
1896 fprintf( io,
"\t%.3e", (*tr).Emis().Aul() * (*tr).EnergyErg() * (*(*tr).Hi()).g() );
1902 else if( (strcmp(chJOB ,
"H2sp" ) == 0) )
1904 fprintf( io,
"PUT SOMETHING HERE!\n" );
1928 return h2.getLine( iElecHi, iVibHi, iRotHi, iElecLo, iVibLo, iRotLo, relint, absint );
1931long int diatomics::getLine(
long iElecHi,
long iVibHi,
long iRotHi,
long iElecLo,
long iVibLo,
long iRotLo,
double *relint,
double *absint )
1941 if( iElecHi!=0 || iElecLo!=0 )
1950 if(
states[ipHi].energy().WN() <
states[ipLo].energy().WN() )
1971 *relint =
H2_SaveLine[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo]/
1980 if(
H2_SaveLine[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo] > 0. )
1982 *absint = log10(
H2_SaveLine[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo]) +
const int FILENAME_PATH_LENGTH_2
NORETURN void BadRead(void)
double RandGauss(double xMean, double s)
double FFmtRead(const char *chCard, long int *ipnt, long int last, bool *lgEOL)
char * read_whole_line(char *chLine, int nChar, FILE *ioIN)
void spsort(realnum x[], long int n, long int iperm[], int kflag, int *ier)
sys_float SDIV(sys_float x)
#define DEBUG_ENTRY(funcname)
realnum & Pelec_esc() const
bool nMatch(const char *chKey) const
double getNumberDefault(const char *chDesc, double fdef)
double getNumberDefaultNegImplLog(const char *chDesc, double fdef)
TransitionProxy::iterator iterator
void H2_PunchDo(FILE *io, char chJOB[], const char chTime[], long int ipPun)
char chH2ColliderLabels[N_X_COLLIDER][chN_X_COLLIDER]
multi_arr< double, 2 > H2_col_rate_out
double Solomon_elec_decay_s
multi_arr< bool, 3 > H2_lgOrtho
double H2_renorm_chemistry
multi_arr< realnum, 3 > H2_X_hminus_formation_distribution
double pops_per_elec[N_ELEC]
long int getLine(long iElecHi, long iVibHi, long iRotHi, long iElecLo, long iVibLo, long iRotLo, double *relint, double *absint)
double Solomon_elec_decay_g
multi_arr< realnum, 3 > H2_dissprob
void set_numLevelsMatrix(long numLevels)
void H2_ReadDissprob(long int nelec)
void H2_Read_hminus_distribution(void)
void H2_ParseSave(Parser &p, char *chHeader)
multi_arr< long int, 2 > ipTransitionSort
valarray< long > ipVib_H2_energy_sort
const double *const dense_total
void H2_Prt_line_tau(void)
multi_arr< realnum, 3 > H2_stat
double Solomon_dissoc_rate_g
double H2_DissocEnergies[N_ELEC]
multi_arr< double, 2 > H2_col_rate_in
multi_arr< realnum, 6 > H2_SaveLine
multi_arr< double, 3 > H2_populations_LTE
valarray< long > ipRot_H2_energy_sort
multi_arr< realnum, 3 > CollRateCoeff
void H2_Prt_column_density(FILE *ioMEAN)
void H2_ReadDissocEnergies(void)
void H2_PrtDepartCoef(void)
void H2_ReadTransprob(long int nelec, TransitionList &trans)
multi_arr< bool, 2 > lgH2_radiative
multi_arr< double, 2 > pops_per_vib
multi_arr< realnum, 2 > H2_X_colden_LTE
void H2_PunchLineStuff(FILE *io, realnum xLimit, long index)
multi_arr< realnum, 2 > H2_X_colden
double Solomon_dissoc_rate_s
valarray< long > nRot_hi[N_ELEC]
multi_arr< double, 3 > H2_old_populations
multi_arr< double, 3 > H2_rad_rate_out
TransitionList::iterator rad_end
multi_arr< realnum, 3 > H2_disske
multi_arr< double, 2 > H2_rad_rate_in
long int nLevels_per_elec[N_ELEC]
void H2_Punch_line_data(FILE *ioPUN, bool lgDoAll)
multi_arr< long int, 3 > ipEnergySort
bool operator<(const level_tmp &second) const
ProxyIterator< qStateConstProxy, qStateConstProxy > const_iterator
ProxyIterator< qStateProxy, qStateConstProxy > iterator
multi_arr< realnum, 3 >::const_iterator mr3ci
FILE * open_data(const char *fname, const char *mode, access_scheme scheme)
realnum GetDopplerWidth(realnum massAMU)
diatomics h2("h2", 4100., &hmi.H2_total, Yan_H2_CS)
void lindst(double xInten, realnum wavelength, const char *chLab, long int ipnt, char chInfo, bool lgOutToo, const char *chComment)
molezone * findspecieslocal(const char buf[])
STATIC char chMolBranch(long iRotHi, long int iRotLo)
long int cdH2_Line(long int iElecHi, long int iVibHi, long int iRotHi, long int iElecLo, long int iVibLo, long int iRotLo, double *relint, double *absint)
static realnum thresh_punline_h2
void prt_wl(FILE *ioOUT, realnum wl)
void prme(const bool lgReset, const TransitionProxy &t)
void Save1Line(const TransitionProxy &t, FILE *io, realnum xLimit, long index, realnum DopplerWidth)
void Save1LineData(const TransitionProxy &t, FILE *io, bool lgCS_2, bool &lgPrint)
t_secondaries secondaries
void PutLine(const TransitionProxy &t, const char *chComment, const char *chLabelTemp)