Home > atmlab > forwardmodel > mpm_iwc93_lwc93_rain89.m

mpm_iwc93_lwc93_rain89

PURPOSE ^

MPM_IWC93_LWC93_RAIN89 Complex refractivity for IWC, LWC and rain.

SYNOPSIS ^

function [nr,alpha] = mpm_iwc93_lwc93_rain89(iwc,lwc,rain,v,t)

DESCRIPTION ^

 MPM_IWC93_LWC93_RAIN89   Complex refractivity for IWC, LWC and rain.

    Returns the complex refractivity following MPM93 for ice water and
    liquid water, and MPM89 for rain. 

    Parameter values for rain rates > 25 mm/h are derived by Christian
    Melsheimer from the Olsen 1978 paper.

 FORMAT   [nr,alpha] = mpm_iwc93_lwc93_rain89(iwc,lwc,rain,v,t)

 OUT      nr     Real part of refractivity [-] (a value just above 1).
          alpha  Absorption [1/m]
 IN       iwc    Ice water content [g/m3]
          lwc    Liquid water content [g/m3]
          rain   Rain rate [mm/h]    
          v      Frequency [Hz]
          t      Temperature [K]

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

DOWNLOAD ^

mpm_iwc93_lwc93_rain89.m

SOURCE CODE ^

0001 % MPM_IWC93_LWC93_RAIN89   Complex refractivity for IWC, LWC and rain.
0002 %
0003 %    Returns the complex refractivity following MPM93 for ice water and
0004 %    liquid water, and MPM89 for rain.
0005 %
0006 %    Parameter values for rain rates > 25 mm/h are derived by Christian
0007 %    Melsheimer from the Olsen 1978 paper.
0008 %
0009 % FORMAT   [nr,alpha] = mpm_iwc93_lwc93_rain89(iwc,lwc,rain,v,t)
0010 %
0011 % OUT      nr     Real part of refractivity [-] (a value just above 1).
0012 %          alpha  Absorption [1/m]
0013 % IN       iwc    Ice water content [g/m3]
0014 %          lwc    Liquid water content [g/m3]
0015 %          rain   Rain rate [mm/h]
0016 %          v      Frequency [Hz]
0017 %          t      Temperature [K]
0018 
0019 % HISTORY: 2003-03-09  Created by Patrick Eriksson
0020 
0021 
0022 function [nr,alpha] = mpm_iwc93_lwc93_rain89(iwc,lwc,rain,v,t)
0023 
0024 
0025 if v < 1e9  |  v > 1000e9
0026   error('This function is only valid for 1 -1000 GHz.');
0027 end
0028 
0029 
0030 %=== Local variables
0031 %
0032 theta = 300 / t;
0033 v_ghz = v/1e9;
0034 
0035 
0036 N = 0;
0037 
0038 
0039 %=== Ice water
0040 %
0041 if iwc < 0
0042   error('Negativ IWC is not allowed.');
0043 elseif iwc > 10
0044   error('IWC > 10 g/m3 is not allowed.');
0045 elseif iwc > 0
0046   eps = eps_ice_liebe93(v,t);
0047   N   = N + iwc*(1.5/0.916)*(eps-1)/(eps+2);
0048 end
0049 
0050 
0051 
0052 %=== Liquid water
0053 %
0054 if lwc < 0
0055   error('Negativ LWC is not allowed.');
0056 elseif lwc > 5
0057   error('LWC > 5 g/m3 is not allowed.');
0058 elseif lwc > 0
0059   eps = eps_water_liebe93(v,t);
0060   N   = N + lwc*1.5*(eps-1)/(eps+2);
0061 end
0062 
0063 
0064 
0065 %=== Rain
0066 %
0067 if rain < 0
0068   error('Negativ rain rate is not allowed.');
0069 elseif rain > 0
0070   %
0071   if rain <= 25
0072     %
0073     if v_ghz <= 2.9
0074       x1 = 3.51e-4;   y1 = 1.03;
0075     elseif v_ghz <= 54
0076       x1 = 2.31e-4;   y1 = 1.42;
0077     elseif v_ghz <= 180
0078       x1 = 0.225;     y1 = -0.301;
0079     else
0080       x1 = 18.6;      y1 = -1.151;
0081     end
0082     %
0083     if v_ghz <= 8.5
0084       x2 = 0.851;     y2 = 0.158;
0085     elseif v_ghz <= 25
0086       x2 = 1.41;      y2 = -0.0779;
0087     elseif v_ghz <= 164
0088       x2 = 2.63;      y2 = -0.272;
0089     else
0090       x2 = 0.616;     y2 = 0.0126;
0091     end
0092     %
0093   else   % rain rate > 25 mm/h
0094     if v > 100e9
0095       error('Frequencies > 100 GHz and rain rates > 25 mm/h are not handled.');
0096     end
0097     %
0098     if v_ghz <= 4.9
0099       x1 = 2.91e-4;   y1 = 0.871;
0100     elseif v_ghz <= 10.7
0101       x1 = 2.76e-5;   y1 = 2.349;
0102     elseif v_ghz <= 40.1
0103       x1 = 1.39e-4;   y1 = 1.668;
0104     elseif v_ghz <= 59.1
0105       x1 = 1.96e-2;   y1 = 0.326;
0106     else
0107       x1 = 0.785;     y1 = -0.578;
0108     end
0109     %
0110     if v_ghz <= 6.2
0111       x2 = 0.911;     y2 = 0.190;
0112     elseif v_ghz <= 23.8
0113       x2 = 1.71;      y2 = -0.156;
0114     elseif v_ghz <= 48.4
0115       x2 = 3.08;      y2 = -0.342;
0116     elseif v_ghz <= 68.2
0117       x2 = 1.28;      y2 = -0.116;
0118     else
0119       x2 = 0.932;     y2 = -4.08e-2;
0120     end
0121     %
0122   end
0123 
0124   fr  = 53 - rain * ( 0.37 - 0.0015 * rain );
0125   y25 = (v_ghz/fr)^2.5;
0126   %
0127   N = N + ...
0128           rain * ( 3.7 - 0.012*rain ) * ( 1- y25/( 1 + y25 ) ) / fr +...
0129           i*x1*v_ghz^y1 * rain^( x2*v_ghz^y2 );
0130 end
0131 
0132 
0133 
0134 %=== Unit conversions
0135 %
0136 nr = 1 + 1e-6 * real( N );  % ppm -> -
0137 %
0138 alpha = 0.182 * v_ghz * imag( N );   % dB/km
0139 alpha = 1e-4 * log(10) * alpha;      % 1/m
0140 
0141 
0142 %alpha = 4 * pi * v * nc/1e6 / constants('SPEED_OF_LIGHT');
0143 
0144

Generated on Mon 15-Sep-2014 13:31:28 by m2html © 2005