.. _program_listing_file_environment_atmosphere.cc:

Program Listing for File atmosphere.cc
======================================

|exhale_lsh| :ref:`Return to documentation for file <file_environment_atmosphere.cc>` (``environment/atmosphere.cc``)

.. |exhale_lsh| unicode:: U+021B0 .. UPWARDS ARROW WITH TIP LEFTWARDS

.. code-block:: cpp

   
   
   #include "lupnt/environment/atmosphere.h"
   
   #include "lupnt/core/constants.h"
   
   namespace lupnt {
   
     double Klobucher(double t_gps, double elevation, double azimuth, double latitude_u,
                      double longitude_u, double freq_Hz, const Vec4d& alpha, const Vec4d& beta) {
       // Compute the earth-centered angle
       double psi = 0.0137 / (elevation + 0.11) - 0.022;
   
       // compute the latitude of the ionospheric pierce point
       double phi_i = latitude_u + psi * std::cos(azimuth);
       if (phi_i > 0.416) {
         phi_i = 0.416;
       } else if (phi_i < -0.416) {
         phi_i = -0.416;
       }
   
       // compute the longitude of the ionospheric pierce point
       double lambda_i = longitude_u + psi * std::sin(azimuth) / std::cos(phi_i);
   
       // compute the geometric latitude of the IPP
       double phi_m = phi_i + 0.064 * std::cos(lambda_i - 1.617);
   
       // compute the local time at the IPP
       double t = 43200 * lambda_i + t_gps;  // sec
   
       // wrap to be within 0 - 86400
       t = fmod(t, 86400);
   
       // compute the amplitude of delay
       double A_I = 0.0;
       for (int i = 0; i < 4; i++) {
         A_I += alpha[i] * std::pow(phi_m, i);
       }
       if (A_I < 0.0) A_I = 0.0;
   
       // compute the period of delay
       double P_I = 0.0;
       for (int i = 0; i < 4; i++) {
         P_I += beta[i] * std::pow(phi_m, i);
       }
       if (P_I < 72000.0) P_I = 72000.0;
   
       // compute the phase of the delay
       double X_I = 2 * PI * (t - 50400) / P_I;
   
       // comute the slant factor
       double F = 1.0 + 16.0 * std::pow(0.53 - elevation, 3);
   
       // compute the ionospheric delay at L1
       double iono_delay = 0.0;
       if (std::abs(X_I) < 1.57) {
         iono_delay = F * (5.0e-9 + A_I * (1 - std::pow(X_I, 2) / 2 + std::pow(X_I, 4) / 24));
       } else {
         iono_delay = F * 5.0e-9;
       }
   
       // convert to the frequency of interest
       iono_delay = iono_delay * std::pow(L1_FREQ / freq_Hz, 2);
       return iono_delay;
     }
   }  // namespace lupnt