OB.DAAC Logo
NASA Logo
Ocean Color Science Software

ocssw V2022
met_cvt.c
Go to the documentation of this file.
1 #include "met_cvt.h"
2 #include <genutils.h>
3 
4 /*
5 
6  met_cvt.c is a collection of meteorological conversion routines
7 
8  Modification history:
9  Programmer Date Description of change
10  ---------- ---- ---------------------
11  W. Robinson 17-Aug-2006 Original development
12  W. Robinson 5 Aug 2013 change to met_cvt
13  */
14 
15 int met_cvt_q_to_rh(int nval, float *pres, int p_type, float *temp,
16  int t_type, float *q, int q_type, float *rh)
17 /*******************************************************************
18 
19  met_cvt_q_to_rh
20 
21  purpose: convert q (specific humidity) into relative humidity
22 
23  Returns type: int - 0 if all is OK
24 
25  Parameters: (in calling order)
26  Type Name I/O Description
27  ---- ---- --- -----------
28  int nval I # of values in following
29  arrays to convert
30  float * pres I pressure array
31  int p_type I pressure type see met_convert.h
32  for the types
33  float * temp I temperature array
34  int t_type I temperature type
35  float * q I specific humidity array
36  int q_type I specific humidity type
37  float * rh O relative humidity in %
38 
39  Modification history:
40  Programmer Date Description of change
41  ---------- ---- ---------------------
42  W. Robinson 12-May-2009 Original development
43 
44  *******************************************************************/ {
45  int ival;
46  double x_h20; /* volume mixing ratio */
47  double es, es_wet, es_ice; /* saturation vapor pressure */
48  double p_lcl, t_lcl, q_lcl, rh_lcl, t_es;
49  /*
50  * loop through the inputs
51  */
52  for (ival = 0; ival < nval; ival++) {
53  /*
54  * get local parameter values and convert them to
55  * hPa (p), C (t), and kg kg^-1 (q)
56  */
57  p_lcl = *(pres + ival);
58  t_lcl = *(temp + ival);
59  q_lcl = *(q + ival);
60  if ((p_lcl != BAD_FLT) && (t_lcl != BAD_FLT) &&
61  (q_lcl != BAD_FLT)) {
62 
63  p_lcl = met_cvt_p_cvt(p_lcl, p_type, MET_UNITS__P_HPA);
64  t_lcl = met_cvt_t_cvt(t_lcl, t_type, MET_UNITS__T_C);
65  q_lcl = met_cvt_q_cvt(q_lcl, q_type, MET_UNITS__Q_KG_KG);
66  /*
67  * then, create the RH
68  * Note that for T < -50 C, the es is not good. To have something
69  * compute es at -50 if t is lower
70  */
71  t_es = (t_lcl < -50.) ? -50. : t_lcl;
72  es_wet = es_coef_wet[0] + t_es * (es_coef_wet[1] + t_es *
73  (es_coef_wet[2] + t_es *
74  (es_coef_wet[3] + t_es * (es_coef_wet[4] + t_es *
75  (es_coef_wet[5] + t_es * es_coef_wet[6])))));
76  es_ice = es_coef_ice[0] + t_es * (es_coef_ice[1] + t_es *
77  (es_coef_ice[2] + t_es *
78  (es_coef_ice[3] + t_es * (es_coef_ice[4] + t_es *
79  (es_coef_ice[5] + t_es * es_coef_ice[6])))));
80  es = (es_wet < es_ice) ? es_wet : es_ice;
81 
82  x_h20 = q_lcl * M_DRY / (q_lcl * M_DRY + (1. - q_lcl) * M_WET);
83 
84  rh_lcl = 100. * x_h20 * p_lcl / es;
85  /*
86  * limit the RH that comes out and transfer
87  */
88  rh_lcl = (rh_lcl < 0) ? 0 : rh_lcl;
89  rh_lcl = (rh_lcl > 100.) ? 100. : rh_lcl;
90 
91  *(rh + ival) = rh_lcl;
92  } else {
93  *(rh + ival) = BAD_FLT;
94  }
95  }
96  return 0;
97 }
98 
99 int met_cvt_rh_to_q( int nval, float *pres, int p_type, float *temp,
100  int t_type, float *rh, float *q, int q_type)
101 /*******************************************************************
102 
103  met_cvt_rh_to_q
104 
105  purpose: convert rh (relative humidity) into specific humidity
106 
107  Returns type: int - 0 if all is OK
108 
109  Parameters: (in calling order)
110  Type Name I/O Description
111  ---- ---- --- -----------
112  int nval I # of values in following
113  arrays to convert
114  float * pres I pressure array
115  int p_type I pressure type see met_convert.h
116  for the types
117  float * temp I temperature array
118  int t_type I temperature type
119  float * rh O relative humidity in %
120  float * q I specific humidity array
121  int q_type I specific humidity type
122 
123  Modification history:
124  Programmer Date Description of change
125  ---------- ---- ---------------------
126  W. Robinson 22-Jan-2019 Original development
127 
128  *******************************************************************/ {
129  int ival;
130  double x_h20; /* volume mixing ratio */
131  double es, es_wet, es_ice; /* saturation vapor pressure */
132  double p_lcl, t_lcl, q_lcl, rh_lcl, t_es, ph2o;
133  /*
134  * loop through the inputs
135  */
136  for (ival = 0; ival < nval; ival++) {
137  /*
138  * get local parameter values and convert them to
139  * hPa (p), C (t), and kg kg^-1 (q)
140  */
141  p_lcl = *(pres + ival);
142  t_lcl = *(temp + ival);
143  rh_lcl = *(rh + ival);
144  if ((p_lcl != BAD_FLT) && (t_lcl != BAD_FLT) &&
145  (rh_lcl != BAD_FLT)) {
146 
147  p_lcl = met_cvt_p_cvt(p_lcl, p_type, MET_UNITS__P_HPA);
148  t_lcl = met_cvt_t_cvt(t_lcl, t_type, MET_UNITS__T_C);
149  rh_lcl = ( rh_lcl > 100. ) ? 100. : rh_lcl;
150  rh_lcl = ( rh_lcl < 0. ) ? 0. : rh_lcl;
151  /*
152  * then, create the Q
153  * Note that for T < -50 C, the es is not good. To have something
154  * compute es at -50 if t is lower
155  */
156  t_es = (t_lcl < -50.) ? -50. : t_lcl;
157  es_wet = es_coef_wet[0] + t_es * (es_coef_wet[1] + t_es *
158  (es_coef_wet[2] + t_es *
159  (es_coef_wet[3] + t_es * (es_coef_wet[4] + t_es *
160  (es_coef_wet[5] + t_es * es_coef_wet[6])))));
161  es_ice = es_coef_ice[0] + t_es * (es_coef_ice[1] + t_es *
162  (es_coef_ice[2] + t_es *
163  (es_coef_ice[3] + t_es * (es_coef_ice[4] + t_es *
164  (es_coef_ice[5] + t_es * es_coef_ice[6])))));
165  es = (es_wet < es_ice) ? es_wet : es_ice;
166 
167  ph2o = es * rh_lcl / 100.;
168 
169  x_h20 = ph2o / p_lcl;
170 
171  q_lcl = x_h20 * M_WET / ( x_h20 * M_WET + ( 1. - x_h20 ) * M_DRY );
172  q_lcl = met_cvt_q_cvt( q_lcl, MET_UNITS__Q_KG_KG, q_type );
173 
174  *(q + ival) = q_lcl;
175  } else {
176  *(q + ival) = BAD_FLT;
177  }
178  }
179  return 0;
180 }
181 
182 int met_cvt_ttd_to_rh(int nval, float *temp, int t_type, float *dwp_temp,
183  int dwp_type, float *rh)
184 /*******************************************************************
185 
186  met_cvt_ttd_to_rh
187 
188  purpose: convert temperature and dewpoint temp into relative humidity
189 
190  Uses the Magnus formula for es from:
191  Reference "The relationship between relative humidity and the dewpoint
192  temperature in moist air", Mark G. Lawrence, BAMS, Feb, 2005,
193  pp. 225 - 233, DOI:10.1175/BAMS-86-2-225
194 
195  Returns type: int - 0 if all is OK
196 
197  Parameters: (in calling order)
198  Type Name I/O Description
199  ---- ---- --- -----------
200  int nval I # of values in following
201  arrays to convert
202  float * temp I temperature array
203  int t_type I temperature type
204  float * dwp_temp I dewpoint temperature
205  int dwp_type I dewpoint type
206  float * rh O relative humidity in %
207 
208  Modification history:
209  Programmer Date Description of change
210  ---------- ---- ---------------------
211  W. Robinson 5 Aug 2013 Original development
212 
213  *******************************************************************/ {
214  int ival;
215  float rh_lcl;
216  double td_lcl, t_lcl;
217  for (ival = 0; ival < nval; ival++) {
218  t_lcl = (double) *(temp + ival);
219  td_lcl = (double) *(dwp_temp + ival);
220  /*
221  * the equations use t in C, so make sure of this
222  */
223  t_lcl = met_cvt_t_cvt(t_lcl, t_type, MET_UNITS__T_C);
224  td_lcl = met_cvt_t_cvt(td_lcl, dwp_type, MET_UNITS__T_C);
225  rh_lcl = 100. * (exp((MAGNUS_A1 * td_lcl) / (MAGNUS_B1 + td_lcl))
226  / exp((MAGNUS_A1 * t_lcl) / (MAGNUS_B1 + t_lcl)));
227 
228  rh_lcl = (rh_lcl < 0) ? 0 : rh_lcl;
229  rh_lcl = (rh_lcl > 100.) ? 100. : rh_lcl;
230 
231  *(rh + ival) = rh_lcl;
232  }
233  return 0;
234 }
235 
236 double met_cvt_p_cvt(double val, int in_type, int out_type)
237 /*******************************************************************
238 
239  met_cvt_p_cvt
240 
241  purpose: convert units for pressure
242 
243  Returns double of new pressure
244 
245  Parameters: (in calling order)
246  Type Name I/O Description
247  ---- ---- --- -----------
248  double val I input pressure
249  int in_type I type of input pressure
250  int out_type I type to convert to
251 
252  Modification history:
253  Programmer Date Description of change
254  ---------- ---- ---------------------
255  W. Robinson 12-May-2009 Original development
256 
257  *******************************************************************/ {
258  float mult;
259 
260  /* assume HPa is normal */
261  if (in_type == out_type)
262  return val;
263  else {
264  switch (in_type) {
265  case MET_UNITS__P_PA:
266  mult = .01;
267  break;
268  case MET_UNITS__P_HPA:
269  mult = 1.;
270  break;
271  }
272  switch (out_type) {
273  case MET_UNITS__P_PA:
274  mult *= 100.;
275  break;
276  case MET_UNITS__P_HPA:
277  mult *= 1.;
278  break;
279  }
280  return mult * val;
281  }
282 }
283 
284 double met_cvt_t_cvt(double val, int in_type, int out_type)
285 /*******************************************************************
286 
287  met_cvt_t_cvt
288 
289  purpose: convert units for temperature
290 
291  Returns double of new temperature
292 
293  Parameters: (in calling order)
294  Type Name I/O Description
295  ---- ---- --- -----------
296  double val I input temperature
297  int in_type I type of input temperature
298  int out_type I type to convert to
299 
300  Modification history:
301  Programmer Date Description of change
302  ---------- ---- ---------------------
303  W. Robinson 12-May-2009 Original development
304 
305  *******************************************************************/ {
306  if (in_type == out_type)
307  return val;
308  else if (in_type == MET_UNITS__T_K)
309  return val - C_IN_K;
310  else
311  return val + C_IN_K;
312 }
313 
314 double met_cvt_q_cvt(double val, int in_type, int out_type)
315 /*******************************************************************
316 
317  met_cvt_q_cvt
318 
319  purpose: convert units for specific humidity
320 
321  Returns double of new specific humidity
322 
323  Parameters: (in calling order)
324  Type Name I/O Description
325  ---- ---- --- -----------
326  double val I input specific humidity
327  int in_type I type of input specific humidity
328  int out_type I type to convert to
329 
330  Modification history:
331  Programmer Date Description of change
332  ---------- ---- ---------------------
333  W. Robinson 12-May-2009 Original development
334 
335  *******************************************************************/ {
336  if (in_type == out_type)
337  return val;
338  else if (in_type == MET_UNITS__Q_KG_KG)
339  return val / 1000.;
340  else
341  return val * 1000.;
342 }
data_t q
Definition: decode_rs.h:74
#define MAGNUS_B1
Definition: met_cvt.h:40
#define MET_UNITS__T_K
Definition: met_cvt.h:26
#define M_DRY
Definition: met_cvt.h:34
#define M_WET
Definition: met_cvt.h:35
#define MET_UNITS__P_PA
Definition: met_cvt.h:24
double met_cvt_t_cvt(double val, int in_type, int out_type)
Definition: met_cvt.c:284
double met_cvt_q_cvt(double val, int in_type, int out_type)
Definition: met_cvt.c:314
integer, parameter double
int met_cvt_q_to_rh(int nval, float *pres, int p_type, float *temp, int t_type, float *q, int q_type, float *rh)
Definition: met_cvt.c:15
#define MET_UNITS__T_C
Definition: met_cvt.h:27
#define BAD_FLT
Definition: jplaeriallib.h:19
int met_cvt_ttd_to_rh(int nval, float *temp, int t_type, float *dwp_temp, int dwp_type, float *rh)
Definition: met_cvt.c:182
#define C_IN_K
Definition: met_cvt.h:36
#define MET_UNITS__P_HPA
Definition: met_cvt.h:25
double met_cvt_p_cvt(double val, int in_type, int out_type)
Definition: met_cvt.c:236
#define MAGNUS_A1
Definition: met_cvt.h:39
int met_cvt_rh_to_q(int nval, float *pres, int p_type, float *temp, int t_type, float *rh, float *q, int q_type)
Definition: met_cvt.c:99
#define MET_UNITS__Q_KG_KG
Definition: met_cvt.h:28
msiBandIdx val
Definition: l1c_msi.cpp:34