一些常用的物理量计算

气象领域中一些常用的物理量计算代码，如水汽压、饱和水汽压、混合比、位温、假相当位温、比湿、饱和比湿、露点温度、水的汽化潜热、相对湿度、绝对湿度等。以及球面距离公式。

---

Formula.java

---

001 /**
002   PACKAGE  : cma.common.atmos
003 
004        功能：气象领域的常用物理量计算
005      文件名：Formula.java
006    相关文件：
007      编写者：刘泽军
008    开始日期：
009    更新日期：
010        备注：
011 
012  */
013 
014 package cma.common.atmos;
015 
016 import java.io.*;
017 import java.lang.*;
018 import java.util.*;
019 
020 public class Formula {
021 /**
022  * 构造函数
023  */
024     public Formula() {
025     }
026     public static double VALUE  = 9999;//Algorithm.DefaultValue;    //缺省值
027 /**
028  * 功能：计算两点距离
029  * 参数：
030  * lon1, lat1   - 第1点位置(度)
031  * lon2, lat2   - 第2点位置(度)
032  * 返回值：
033  *      两点距离
034  */
035     public static double distanceDegree(double lon1, double lat1, double lon2, double lat2) {
036         return(distance(Math.toRadians(lon1), Math.toRadians(lat1), Math.toRadians(lon1), Math.toRadians(lat1)));
037     }
038 /**
039  * 功能：计算两点距离
040  * 参数：
041  * lon1, lat1   - 第1点位置(弧度)
042  * lon2, lat2   - 第2点位置(弧度)
043  * 返回值：
044  *      两点距离
045  */
046     public static double distance(double lon1, double lat1, double lon2, double lat2) {//(弧度)
047         return(Meteorology.RADIUS * Math.acos(Math.sin(lat1)*Math.sin(lat2) + Math.cos(lat1)*Math.cos(lat2)*Math.cos(lon2-lon1)));
048         //地球半径*arccsos(sin(wd1)sin(wd2)+cos(wd1)cos(wd2)cos(jd2-jd1)) 球面距离？
049     }
050     public static double w(double _p/*气压*/, double _e/*水汽压*/) {//混合比
051         if( VALUE != _p && VALUE != _e && 0 != _p && _p != _e ) {
052             return(0.62198 * _e / ( _p - _e ));
053         }
054         return(VALUE);
055     }
056     public static double e(double _p/*气压*/, double _q/*比湿*/) {//水汽压
057         if( VALUE != _p && VALUE != _q && 0 != _p ) {
058             return(_p * _q / 0.62198);
059         }
060         return(VALUE);
061     }
062     public static double e(double _td/*露点温度*/) {//水汽压
063         //《大气科学常用公式》P382 td => e
064         if( _td == 0.0 ) {
065             return(6.1078);
066         }
067         else if( VALUE != _td && -237.3 != _td ) {
068             return(6.1078 * Math.pow(10, 7.5/(237.3/_td+1)));
069         }
070         else {
071             return(VALUE);
072         }
073     }
074     public static double E(double _t/*温度*/) {//饱和水汽压
075         //《大气科学常用公式》P2
076         if( VALUE != _t && -243.92 != _t ) {
077             return(6.1078 * Math.pow(10,7.69*_t/(_t+243.92)));
078         }
079         return(VALUE);
080     }
081     public static double a(double _e/*水汽压*/, double _t/*温度*/) {//绝对湿度=e/Rv/T
082         double  _Rv  = 461.51;//J . kg-1 . K-1
083         if( VALUE != _e && VALUE != _t && -273.16 != _t ) {
084             return(_e / (_Rv*(273.16+_t)));
085         }
086         return(VALUE);
087     }
088     public static double r(double _e/*水汽压*/, double _E/*饱和水汽压*/) {//相对湿度=e/E*100%
089         if( VALUE != _e && VALUE != _E ) {
090             return(_e / _E * 100.0);
091         }
092         return(VALUE);
093     }
094     public static double theta(double _p/*气压*/, double _t/*温度*/) {//理想气体的位温
095         //以干空气的位温代替
096         return(thetad(_p, _t));
097     }
098     public static double thetad(double _p/*气压*/, double _t/*温度*/) {//干空气的位温
099         if( VALUE != _t && VALUE != _p && 0.0 != _p ) {
100             return( (_t + 273.16) * Math.pow(1000/_p, 0.286) );//《大气科学常用公式》P466 pd= p
101         }
102         return(VALUE);
103     }
104     public static double thetam(double _p/*气压*/, double _t/*温度*/, double _e/*水汽压*/, double _q/*比湿*/) {//湿空气的位温
105         double  _Rm      = 287.04 * ( 1 + 0.608 * _q );
106         double  _Cpm     = 1005 * ( 1 + 0.86 * _q );
107         if( VALUE != _t && VALUE != _p && 0.0 != _p && VALUE != _e && 0.0 != _p + _e && VALUE != _q ) {
108             return( (_t + 273.16) * Math.pow(1000/(_p+_e),_Rm/_Cpm) );//《大气科学常用公式》P466
109         }
110         return(VALUE);
111     }
112     public static double blowupTheta1(double _p/*气压*/, double _t/*温度*/, double _td/*露点温度*/) {//溃变 - θ1
113         return(theta(_p, _t));
114     }
115     public static double blowupTheta2(double _p/*气压*/, double _t/*温度*/, double _td/*露点温度*/) {//溃变 - θ2
116         if( _p == 0 || _p == VALUE || _t == VALUE || _td == VALUE )  return(VALUE);
117         double  _o1 = blowupTheta1(_p, _t, _td);
118         if( VALUE == _o1 ) {
119             return(VALUE);
120         }
121         double  _h  = 124 * ( _t - _td );
122         double  _tc = _t - 0.0098 * _h;
123         double  _e  = 6.1078 * Math.pow( 10, (7.5*_tc)/(237.3+_tc) );
124         double  _qs = 0.62198 * _e / _p;
125         double  _cp = 1005;
126         double  _l  = 2500000;
127         double  _o2 = _o1 * Math.exp( (_l*_qs) / (_cp*(_tc+273.16)) );
128         return(_o2);
129     }
130     public static double  blowupTheta3(double _p/*气压*/, double _t/*温度*/, double _td/*露点温度*/) {//溃变 - θ3
131         if( _p == 0 || _p == VALUE || _t == VALUE || _td == VALUE )  return(VALUE);
132         double  _o1 = blowupTheta1(_p, _t, _td);
133         if( _o1 == VALUE ) return(VALUE);
134         double  _e1 = 6.1078 * Math.pow( 10, (7.5*_t)/(237.3+_t) );
135         double  _cp = 1005;
136         double  _q  = 0.62198 * _e1 / _p;
137         double  _l  = 2500000;
138         double  _o3 = _o1 * Math.exp((_l*_q)/(_cp*(_t+273.16)));
139         return( _o3 );
140     }
141     public static double  thetase(double _p/*气压*/, double _t/*温度*/) {//假相当位温
142         double  thetase = VALUE;
143         return(thetase);
144     }
145     public static double  q(double _e/*水汽压*/, double _p/*气压*/) {//比湿(e为饱和水汽压返回饱和比湿，否则返回比湿)
146         if( VALUE != _e && VALUE != _p && 0.0 != _p ) {
147             return(0.62198 * _e / _p);//《大气科学常用公式》P382
148         }
149         return(VALUE);
150     }
151     public static double qs(double _E/*饱和水汽压*/, double _p/*气压*/) {//饱和比湿
152         return(q(_E, _p));
153     }
154     public static double dv(double _dd/*风向*/, double _ff/*风速*/) {//风的垂直分量
155         double  vv = VALUE;                                             //    ^
156         if( VALUE != _dd && VALUE != _ff ) {                            //   /|/N
157                                                                         //W   |   E
158             vv = Math.abs( _ff * Math.cos(Math.toRadians(_dd)));        //----+--->
159             if( _dd >= 270.0 || _dd < 90.0 ) {                          //    |
160                 vv = -1.0 * Math.abs(vv);                               //    | S
161             }                                                           //    |
162         }
163         return(vv);
164     }
165     public static double du(double _dd/*风向*/, double _ff/*风速*/) {//风的水平分量
166         double  _vh = VALUE;                                            //    ^
167         if( VALUE != _dd && VALUE != _ff ) {                            //   /|/N
168                                                                         //W   |   E
169             _vh = Math.abs( _ff * Math.sin(Math.toRadians(_dd)) );      //----+--->
170             if( _dd <= 180.0 ) {                                        //    |
171                 _vh = - Math.abs( _vh );                                //    | S
172             }                                                           //    |
173         }
174         return(_vh);
175     }
176     public static double  td(double _e/*水汽压*/) {//利用水汽压计算露点温度
177         double  _td = VALUE;
178         double  _A  = 6.1078, _a = 7.5, _b = 237.3;
179         if( _e == VALUE ) {
180             _td = VALUE;
181         }
182         else if( _e == _A ) {
183             _td = 0;
184         }
185         else {
186             _td = _b / ( _a / Math.log(_e/_A)/Math.log(10.0) - 1);//《大气科学常用公式》P382
187         }
188         return(_td);
189     }
190     public static double  td(double _p/*气压*/,double _t/*温度*/,double _r/*相对湿度*/) {//利用气压、气温及相对湿度计算露点温度 - 待检查
191         double  _E   = E(_t);// t => E
192         double  _e   = _E == VALUE ? VALUE : _E * _r / 100.0;// 《大气科学常用公式》P381 r = e / E * 100 %
193         return( td(_e) );
194     }
195     public static double  lv(double _t/*温度*/)//水的汽化潜热
196 {//单位：J . kg -1
197         double  _Rv  = 461.51;//J . kg-1 . K -1
198         double  _a   = 7.726,    _a1  = 7.44;
199         double  _b   = 244.88,   _b1  = 235.33;
200         double  _lv  = VALUE;
201         if( _t != VALUE && _t <= 0 && _t != -_b ) {
202             _lv = _a * _b * _Rv * Math.log(10) * Math.pow( (_t + 237.16)/(_b+_t), 2 );
203         }
204         else if( _t != VALUE && _t > 0 && _t != -_b1 ) {
205             _lv = _a1 * _b1 * _Rv * Math.log(10) * Math.pow( (_t + 237.16)/(_b1+_t), 2 );
206         }
207         return(_lv);
208     }
209 /*
210     public static String LongLat2dfm(double value) {
211         double  d1,f1,m1;
212         TStringList *StringList1 = new TStringList();
213         TStringList *StringList2 = new TStringList();
214         FountainExtractString(FormatFloat("0.000000",value),".",false,StringList1);
215         FountainExtractString(FormatFloat("0.000000",(AnsiString("0.")+StringList1->Strings[1]).ToDouble()*60),".",false,StringList2);
216 
217         d1 = StringList1->Strings[0].ToDouble();
218         f1 = StringList2->Strings[0].ToDouble();
219         m1 = (AnsiString("0.")+StringList2->Strings[1]).ToDouble()*60;
220 
221         if( m1 >= 59.5 ) {
222             m1 = 0.0;
223             f1 += 1.0;
224         }
225         if( f1 >= 59.5 ) {
226             f1 = 0.0;
227             d1 += 1.0;
228         }
229 
230         AnsiString FReturnStr =
231                 FormatFloat("0",d1)+"°"+
232                 FormatFloat("00",f1)+"′"+
233                 FormatFloat("00",m1)+"″";
234         delete StringList1;
235         delete StringList2;
236         return(FReturnStr);
237     }
238 */
239 }