Changeset 1710


Ignore:
Timestamp:
03/17/10 15:37:26 (14 years ago)
Author:
Max Voronkov
Message:

added the code to calculate dry refractivity of the atmosphere

Location:
trunk/src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/STAtmosphere.cpp

    r1709 r1710  
    173173}
    174174
     175/**
     176 * Compute the complex refractivity of the dry components of the atmosphere
     177 * (oxygen lines) at the given frequency.
     178 * @param[in] freq frequency (Hz)
     179 * @param[in] temperature air temperature (K)
     180 * @param[in] pDry partial pressure of dry components (Pascals)
     181 * @param[in] pVapour partial pressure of water vapour (Pascals)
     182 * @return complex refractivity
     183 **/
     184std::complex<double> STAtmosphere::dryRefractivity(double freq, double temperature,
     185                     double pDry, double pVapour)
     186{
     187  // the number of parameters per atmospheric line and the number of lines taken into account
     188  const size_t nLineParams = 7;
     189  const size_t nLines = 48;
     190  // actual tabulated values
     191  const double lines[nLines][nLineParams] =
     192    {{49.452379,    0.12E-6, 11.830,  8.40E-3, 0.0,  5.60E-3,  1.7},
     193     {49.962257,    0.34E-6, 10.720,  8.50E-3, 0.0,  5.60E-3,  1.7},
     194     {50.474238,    0.94E-6,  9.690,  8.60E-3, 0.0,  5.60E-3,  1.7},
     195     {50.987748,    2.46E-6,  8.690,  8.70E-3, 0.0,  5.50E-3,  1.7},
     196     {51.503350,    6.08E-6,  7.740,  8.90E-3, 0.0,  5.60E-3,  1.8},
     197     {52.021409,   14.14E-6,  6.840,  9.20E-3, 0.0,  5.50E-3,  1.8},
     198     {52.542393,   31.02E-6,  6.000,  9.40E-3, 0.0,  5.70E-3,  1.8},
     199     {53.066906,   64.10E-6,  5.220,  9.70E-3, 0.0,  5.30E-3,  1.9},
     200     {53.595748,  124.70E-6,  4.480, 10.00E-3, 0.0,  5.40E-3,  1.8},
     201     {54.129999,  228.00E-6,  3.810, 10.20E-3, 0.0,  4.80E-3,  2.0},
     202     {54.671157,  391.80E-6,  3.190, 10.50E-3, 0.0,  4.80E-3,  1.9},
     203     {55.221365,  631.60E-6,  2.620, 10.79E-3, 0.0,  4.17E-3,  2.1},
     204     {55.783800,  953.50E-6,  2.115, 11.10E-3, 0.0,  3.75E-3,  2.1},
     205     {56.264777,  548.90E-6,  0.010, 16.46E-3, 0.0,  7.74E-3,  0.9},
     206     {56.363387, 1344.00E-6,  1.655, 11.44E-3, 0.0,  2.97E-3,  2.3},
     207     {56.968180, 1763.00E-6,  1.255, 11.81E-3, 0.0,  2.12E-3,  2.5},
     208     {57.612481, 2141.00E-6,  0.910, 12.21E-3, 0.0,  0.94E-3,  3.7},
     209     {58.323874, 2386.00E-6,  0.621, 12.66E-3, 0.0, -0.55E-3, -3.1},
     210     {58.446589, 1457.00E-6,  0.079, 14.49E-3, 0.0,  5.97E-3,  0.8},
     211     {59.164204, 2404.00E-6,  0.386, 13.19E-3, 0.0, -2.44E-3,  0.1},
     212     {59.590982, 2112.00E-6,  0.207, 13.60E-3, 0.0,  3.44E-3,  0.5},
     213     {60.306057, 2124.00E-6,  0.207, 13.82E-3, 0.0, -4.13E-3,  0.7},
     214     {60.434775, 2461.00E-6,  0.386, 12.97E-3, 0.0,  1.32E-3, -1.0},
     215     {61.150558, 2504.00E-6,  0.621, 12.48E-3, 0.0, -0.36E-3,  5.8},
     216     {61.800152, 2298.00E-6,  0.910, 12.07E-3, 0.0, -1.59E-3,  2.9},
     217     {62.411212, 1933.00E-6,  1.255, 11.71E-3, 0.0, -2.66E-3,  2.3},
     218     {62.486253, 1517.00E-6,  0.078, 14.68E-3, 0.0, -4.77E-3,  0.9},
     219     {62.997974, 1503.00E-6,  1.660, 11.39E-3, 0.0, -3.34E-3,  2.2},
     220     {63.568515, 1087.00E-6,  2.110, 11.08E-3, 0.0, -4.17E-3,  2.0},
     221     {64.127764,  733.50E-6,  2.620, 10.78E-3, 0.0, -4.48E-3,  2.0},
     222     {64.678900,  463.50E-6,  3.190, 10.50E-3, 0.0, -5.10E-3,  1.8},
     223     {65.224067,  274.80E-6,  3.810, 10.20E-3, 0.0, -5.10E-3,  1.9},
     224     {65.764769,  153.00E-6,  4.480, 10.00E-3, 0.0, -5.70E-3,  1.8},
     225     {66.302088,   80.09E-6,  5.220,  9.70E-3, 0.0, -5.50E-3,  1.8},
     226     {66.836827,   39.46E-6,  6.000,  9.40E-3, 0.0, -5.90E-3,  1.7},
     227     {67.369595,   18.32E-6,  6.840,  9.20E-3, 0.0, -5.60E-3,  1.8},
     228     {67.900862,    8.01E-6,  7.740,  8.90E-3, 0.0, -5.80E-3,  1.7},
     229     {68.431001,    3.30E-6,  8.690,  8.70E-3, 0.0, -5.70E-3,  1.7},
     230     {68.960306,    1.28E-6,  9.690,  8.60E-3, 0.0, -5.60E-3,  1.7},
     231     {69.489021,    0.47E-6, 10.720,  8.50E-3, 0.0, -5.60E-3,  1.7},
     232     {70.017342,    0.16E-6, 11.830,  8.40E-3, 0.0, -5.60E-3,  1.7},
     233     {118.750341,  945.00E-6,  0.000, 15.92E-3, 0.0, -0.44E-3,  0.9},
     234     {368.498350,   67.90E-6,  0.020, 19.20E-3, 0.6,  0.00E00,  1.0},
     235     {424.763120,  638.00E-6,  0.011, 19.16E-3, 0.6,  0.00E00,  1.0},
     236     {487.249370,  235.00E-6,  0.011, 19.20E-3, 0.6,  0.00E00,  1.0},
     237     {715.393150,   99.60E-6,  0.089, 18.10E-3, 0.6,  0.00E00,  1.0},
     238     {773.838730,  671.00E-6,  0.079, 18.10E-3, 0.6,  0.00E00,  1.0},
     239     {834.145330,  180.00E-6,  0.079, 18.10E-3, 0.6,  0.00E00,  1.0}};
     240     
     241  // convert to the units of Liebe
     242  const double theta = 300./temperature;
     243  const double kPaPVap = 0.001*pVapour;
     244  const double kPaPDry = 0.001*pDry;
     245  const double fGHz = freq * 1e-9;
     246 
     247  // some coefficients
     248  const double ap = 1.4e-10*(1-1.2e-5*std::pow(fGHz,1.5));
     249  const double gamma0 = 5.6e-3*(kPaPDry + 1.1*kPaPVap)*std::pow(theta,0.8);
     250  // initial refractivity
     251  std::complex<double> result(2.588*kPaPDry*theta +
     252         3.07e-4*(1.0/(1.0+std::pow(fGHz/gamma0,2))-1)*kPaPDry*theta*theta,
     253         (2*3.07e-4/(gamma0*(1+std::pow(fGHz/gamma0,2))*(1+std::pow(fGHz/60,2))) +
     254          ap*kPaPDry*std::pow(theta,2.5))*fGHz*kPaPDry*theta*theta);
     255  // sum the contributions of all the lines
     256  for (size_t l = 0; l < nLines; ++l) {
     257       const double S = lines[l][1]*kPaPDry*std::pow(theta,3)*exp(lines[l][2]*(1.-theta));
     258       const double gamma = lines[l][3]*(kPaPDry*std::pow(theta,0.8-lines[l][4]) + 1.1*kPaPVap*theta);
     259       const double delta = lines[l][5]*kPaPDry*std::pow(theta,lines[l][6]);
     260       const double x = (lines[l][0]-fGHz)*(lines[l][0]-fGHz) + gamma*gamma;
     261       const double y = (lines[l][0]+fGHz)*(lines[l][0]+fGHz) + gamma*gamma;
     262       const double z = (lines[l][0]+gamma*gamma/lines[l][0]);
     263       result += std::complex<double> (S*( (z-fGHz)/x + (z+fGHz)/y - 2./lines[l][0] +
     264                                  delta*(1/x-1/y)*gamma*fGHz/lines[l][0]),
     265               S*( (1/x+1/y)*gamma*fGHz/lines[l][0] -
     266               delta*((lines[l][0]-fGHz)/x + (lines[l][0]+fGHz)/y)*fGHz/lines[l][0]));       
     267  }
     268 
     269  return result;
     270}
  • trunk/src/STAtmosphere.h

    r1709 r1710  
    4343#define STATMOSPHERE_H
    4444
     45// std includes
    4546#include <vector>
     47#include <complex>
    4648
    4749namespace asap {
     
    133135  static double wvSaturationPressure(double temperature);
    134136   
     137  /**
     138   * Compute the complex refractivity of the dry components of the atmosphere
     139   * (oxygen lines) at the given frequency.
     140   * @param[in] freq frequency (Hz)
     141   * @param[in] temperature air temperature (K)
     142   * @param[in] pDry partial pressure of dry components (Pascals)
     143   * @param[in] pVapour partial pressure of water vapour (Pascals)
     144   * @return complex refractivity
     145   **/
     146   static std::complex<double> dryRefractivity(double freq, double temperature,
     147                     double pDry, double pVapour);
     148   
    135149private:
    136150 
Note: See TracChangeset for help on using the changeset viewer.