00001
00002
00003 #include "cddefines.h"
00004 #include "atmdat.h"
00005 #include "thirdparty.h"
00006 t_atmdat atmdat;
00007
00008 double ****HS_He1_Xsectn;
00009 double ****HS_He1_Energy;
00010 double *****OP_Helike_Xsectn;
00011 double *****OP_Helike_Energy;
00012 long ****OP_Helike_NumPts;
00013
00014 void ReadCollisionRateTable( CollRateCoeffArray& coll_rate_table, FILE* io, FunctPtr GetIndices, long nMolLevs, long nTemps, long nTrans )
00015 {
00016 DEBUG_ENTRY( "ReadCollisionRateTable()" );
00017
00018 char chLine[INPUT_LINE_LENGTH] = "";
00019
00020
00021
00022 ASSERT( nTemps != 0 && nTrans != 0 );
00023
00024
00025 while( read_whole_line( chLine, (int)sizeof(chLine), io ) != NULL )
00026 {
00027 if( chLine[0] == '!' || chLine[0] == '#' || chLine[0] == '\n' )
00028 continue;
00029 else
00030 break;
00031 }
00032 ASSERT( strlen( chLine ) > 0 );
00033
00034
00035 {
00036 char *chColltemp = strtok(chLine," \t\n");
00037 while( chColltemp != NULL )
00038 {
00039 coll_rate_table.temps.push_back( atof(chColltemp) );
00040 chColltemp = strtok(NULL," \t\n");
00041 }
00042
00043 if( nTemps < 0 )
00044 nTemps = coll_rate_table.temps.size();
00045 else
00046 ASSERT( (unsigned)nTemps == coll_rate_table.temps.size() );
00047 }
00048
00049 ASSERT( coll_rate_table.collrates.size() == 0 );
00050 coll_rate_table.collrates.reserve( nMolLevs );
00051 for( long ipHi=0; ipHi<nMolLevs; ipHi++ )
00052 {
00053 coll_rate_table.collrates.reserve( ipHi, nMolLevs );
00054 for( long ipLo=0; ipLo<nMolLevs; ipLo++ )
00055 {
00056 coll_rate_table.collrates.reserve( ipHi, ipLo, nTemps );
00057 }
00058 }
00059 coll_rate_table.collrates.alloc();
00060
00061 coll_rate_table.collrates.zero();
00062
00063 long ipTrans = 0;
00064 while( read_whole_line( chLine, (int)sizeof(chLine), io ) != NULL )
00065 {
00066 if( chLine[0] == '!' || chLine[0] == '#' || chLine[0] == '\n' )
00067 continue;
00068
00069 long i = 1;
00070 long ipHi = -1, ipLo = -1;
00071 (*GetIndices)( ipHi, ipLo, chLine, i );
00072 ipTrans++;
00073
00074
00075 if( ipHi == -1 && ipLo == -1 )
00076 continue;
00077
00078
00079 if( ipLo >= nMolLevs || ipHi >= nMolLevs )
00080 {
00081 if( nTrans > 0 && ipTrans == nTrans )
00082 break;
00083 else
00084 continue;
00085 }
00086
00087
00088 if( ipHi < ipLo )
00089 {
00090 ASSERT( ipLo == nMolLevs - 1);
00091 long temp = ipHi;
00092 ipHi = ipLo;
00093 ipLo = temp;
00094 }
00095
00096 ASSERT( ipHi >= 0 );
00097 ASSERT( ipLo >= 0 );
00098
00099 bool lgEOL = false;
00100 for( long j = 0; j < nTemps; ++j )
00101 {
00102 coll_rate_table.collrates[ipHi][ipLo][j] =
00103 (double)FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );
00104 ASSERT( !lgEOL );
00105 }
00106
00107
00108 FFmtRead( chLine, &i, sizeof(chLine), &lgEOL );
00109 ASSERT( lgEOL );
00110
00111 {
00112 enum {DEBUG_LOC=false};
00113 if( DEBUG_LOC )
00114 {
00115 printf("The values of up and lo are %ld & %ld \n",ipHi,ipLo);
00116 printf("The collision rates are");
00117 for( long i = 0; i < nTemps; ++i )
00118 {
00119 printf( "\n %e", coll_rate_table.collrates[ipHi][ipLo][i]);
00120 }
00121 printf("\n");
00122 }
00123 }
00124
00125 if( nTrans > 0 && ipTrans == nTrans )
00126 break;
00127 }
00128 ASSERT( ipTrans > 0 );
00129 if( nTrans > 0 )
00130 ASSERT( ipTrans == nTrans );
00131
00132 return;
00133 }
00134
00135 double InterpCollRate( const CollRateCoeffArray& rate_table, const long& ipHi, const long& ipLo, const double& ftemp)
00136 {
00137 DEBUG_ENTRY("InterpCollRate()");
00138 double ret_collrate = 0.;
00139
00140 if( rate_table.temps.size() == 0 )
00141 {
00142 return 0.;
00143 }
00144
00145 if( ftemp < rate_table.temps[0] )
00146 {
00147 ret_collrate = rate_table.collrates[ipHi][ipLo][0];
00148 }
00149 else if( ftemp > rate_table.temps.back() )
00150 {
00151 ret_collrate = rate_table.collrates[ipHi][ipLo][ rate_table.temps.size()-1 ];
00152 }
00153 else if( rate_table.temps.size() == 1 )
00154 {
00155
00156 ret_collrate = rate_table.collrates[ipHi][ipLo][0];
00157 }
00158 else
00159 {
00160 ret_collrate = linint(&rate_table.temps[0],
00161 &rate_table.collrates[ipHi][ipLo][0],
00162 rate_table.temps.size(),
00163 ftemp);
00164 }
00165
00166 ASSERT( !isnan( ret_collrate ) );
00167 return(ret_collrate);
00168 }
00169