MPTRAC
Macros | Functions
met_map.c File Reference

Extract map from meteorological data. More...

#include "mptrac.h"

Go to the source code of this file.

Macros

#define NX   EX
 Maximum number of longitudes. More...
 
#define NY   EY
 Maximum number of latitudes. More...
 

Functions

void usage (void)
 Print command-line help. More...
 
int main (int argc, char *argv[])
 

Detailed Description

Extract map from meteorological data.

Definition in file met_map.c.

Macro Definition Documentation

◆ NX

#define NX   EX

Maximum number of longitudes.

Definition at line 32 of file met_map.c.

◆ NY

#define NY   EY

Maximum number of latitudes.

Definition at line 35 of file met_map.c.

Function Documentation

◆ usage()

void usage ( void  )

Print command-line help.

Definition at line 460 of file met_map.c.

461 {
462
463 printf("\nMPTRAC met_map tool.\n\n");
464 printf("Extract maps from meteorological data.\n");
465 printf("\n");
466 printf("Usage:\n");
467 printf(" met_map <ctl> <map.tab> <met0> [<met1> ...]\n");
468 printf("\n");
469 printf("Arguments:\n");
470 printf(" <ctl> Control file.\n");
471 printf(" <map.tab> Output table.\n");
472 printf(" <met*> Meteorological input files.\n");
473 printf("\nFurther information:\n");
474 printf(" Manual: https://slcs-jsc.github.io/mptrac/\n");
475}

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 49 of file met_map.c.

51 {
52
53 ctl_t ctl;
54
55 clim_t *clim;
56
57 met_t *met;
58
59 dd_t *dd;
60
61 FILE *out;
62
63 static double *timem, ps, *psm, ts, *tsm, zs, *zsm, us, *usm, vs, *vsm,
64 ess, *essm, nss, *nssm, shf, *shfm, lsm, *lsmm, sst, *sstm, pbl, *pblm,
65 pt, *ptm, t, *pm, *tm, u, *um, v, *vm, w, *wm, h2o, *h2om, h2ot, *h2otm,
66 o3, *o3m, *hno3m, *ohm, *h2o2m, *ho2m, *o1dm, *tdewm, *ticem, *tnatm,
67 lwc, *lwcm, rwc, *rwcm, iwc, *iwcm, swc, *swcm, cc, *ccm, z, *zm,
68 pv, *pvm, zt, *ztm, tt, *ttm, pct, *pctm, pcb, *pcbm, cl, *clm,
69 plcl, *plclm, plfc, *plfcm, pel, *pelm, cape, *capem, cin, *cinm,
70 o3c, *o3cm, *rhm, *rhicem, ptop, pbot, t0, lon, lons[NX], lat, lats[NY];
71
72 static int *np, *npc, *npt, nx, ny;
73
74 /* Allocate... */
75 ALLOC(clim, clim_t, 1);
76 ALLOC(met, met_t, 1);
77 ALLOC(dd, dd_t, 1);
78 ALLOC(timem, double,
79 NX * NY);
80 ALLOC(psm, double,
81 NX * NY);
82 ALLOC(tsm, double,
83 NX * NY);
84 ALLOC(zsm, double,
85 NX * NY);
86 ALLOC(usm, double,
87 NX * NY);
88 ALLOC(vsm, double,
89 NX * NY);
90 ALLOC(essm, double,
91 NX * NY);
92 ALLOC(nssm, double,
93 NX * NY);
94 ALLOC(shfm, double,
95 NX * NY);
96 ALLOC(lsmm, double,
97 NX * NY);
98 ALLOC(sstm, double,
99 NX * NY);
100 ALLOC(pblm, double,
101 NX * NY);
102 ALLOC(ptm, double,
103 NX * NY);
104 ALLOC(pm, double,
105 NX * NY);
106 ALLOC(tm, double,
107 NX * NY);
108 ALLOC(um, double,
109 NX * NY);
110 ALLOC(vm, double,
111 NX * NY);
112 ALLOC(wm, double,
113 NX * NY);
114 ALLOC(h2om, double,
115 NX * NY);
116 ALLOC(h2otm, double,
117 NX * NY);
118 ALLOC(o3m, double,
119 NX * NY);
120 ALLOC(hno3m, double,
121 NX * NY);
122 ALLOC(ohm, double,
123 NX * NY);
124 ALLOC(h2o2m, double,
125 NX * NY);
126 ALLOC(ho2m, double,
127 NX * NY);
128 ALLOC(o1dm, double,
129 NX * NY);
130 ALLOC(tdewm, double,
131 NX * NY);
132 ALLOC(ticem, double,
133 NX * NY);
134 ALLOC(tnatm, double,
135 NX * NY);
136 ALLOC(lwcm, double,
137 NX * NY);
138 ALLOC(rwcm, double,
139 NX * NY);
140 ALLOC(iwcm, double,
141 NX * NY);
142 ALLOC(swcm, double,
143 NX * NY);
144 ALLOC(ccm, double,
145 NX * NY);
146 ALLOC(zm, double,
147 NX * NY);
148 ALLOC(pvm, double,
149 NX * NY);
150 ALLOC(ztm, double,
151 NX * NY);
152 ALLOC(ttm, double,
153 NX * NY);
154 ALLOC(pctm, double,
155 NX * NY);
156 ALLOC(pcbm, double,
157 NX * NY);
158 ALLOC(clm, double,
159 NX * NY);
160 ALLOC(plclm, double,
161 NX * NY);
162 ALLOC(plfcm, double,
163 NX * NY);
164 ALLOC(pelm, double,
165 NX * NY);
166 ALLOC(capem, double,
167 NX * NY);
168 ALLOC(cinm, double,
169 NX * NY);
170 ALLOC(o3cm, double,
171 NX * NY);
172 ALLOC(rhm, double,
173 NX * NY);
174 ALLOC(rhicem, double,
175 NX * NY);
176 ALLOC(np, int,
177 NX * NY);
178 ALLOC(npc, int,
179 NX * NY);
180 ALLOC(npt, int,
181 NX * NY);
182
183 /* Print usage information... */
184 USAGE;
185
186 /* Check arguments... */
187 if (argc < 4)
188 ERRMSG("Missing or invalid command-line arguments.\n\n"
189 "Usage: met_map <ctl> <map.tab> <met0> [<met1> ...]\n\n"
190 "Use -h for full help.");
191
192 /* Read control parameters... */
193 mptrac_read_ctl(argv[1], argc, argv, &ctl);
194 double p0 = P(scan_ctl(argv[1], argc, argv, "MAP_Z0", -1, "10", NULL));
195 double lon0 = scan_ctl(argv[1], argc, argv, "MAP_LON0", -1, "-180", NULL);
196 double lon1 = scan_ctl(argv[1], argc, argv, "MAP_LON1", -1, "180", NULL);
197 double dlon = scan_ctl(argv[1], argc, argv, "MAP_DLON", -1, "-999", NULL);
198 double lat0 = scan_ctl(argv[1], argc, argv, "MAP_LAT0", -1, "-90", NULL);
199 double lat1 = scan_ctl(argv[1], argc, argv, "MAP_LAT1", -1, "90", NULL);
200 double dlat = scan_ctl(argv[1], argc, argv, "MAP_DLAT", -1, "-999", NULL);
201 double theta = scan_ctl(argv[1], argc, argv, "MAP_THETA", -1, "-999", NULL);
202
203 /* Read climatological data... */
204 mptrac_read_clim(&ctl, clim);
205
206 /* Loop over files... */
207 for (int i = 3; i < argc; i++) {
208
209 /* Read meteorological data... */
210 if (!mptrac_read_met(argv[i], &ctl, clim, met, dd))
211 continue;
212
213 /* Set horizontal grid... */
214 if (dlon <= 0)
215 dlon = fabs(met->lon[1] - met->lon[0]);
216 if (dlat <= 0)
217 dlat = fabs(met->lat[1] - met->lat[0]);
218 if (lon0 < -360 && lon1 > 360) {
219 lon0 = gsl_stats_min(met->lon, 1, (size_t) met->nx);
220 lon1 = gsl_stats_max(met->lon, 1, (size_t) met->nx);
221 }
222 nx = ny = 0;
223 for (lon = lon0; lon <= lon1 + 0.001; lon += dlon) {
224 lons[nx] = round(lon * 1e3) / 1e3;
225 if ((++nx) >= NX)
226 ERRMSG("Too many longitudes!");
227 }
228 if (lat0 < -90 && lat1 > 90) {
229 lat0 = gsl_stats_min(met->lat, 1, (size_t) met->ny);
230 lat1 = gsl_stats_max(met->lat, 1, (size_t) met->ny);
231 }
232 for (lat = lat0; lat <= lat1 + 0.001; lat += dlat) {
233 lats[ny] = round(lat * 1e3) / 1e3;
234 if ((++ny) >= NY)
235 ERRMSG("Too many latitudes!");
236 }
237
238 /* Average... */
239 for (int ix = 0; ix < nx; ix++)
240 for (int iy = 0; iy < ny; iy++) {
241
242 /* Find pressure level for given theta level... */
243 INTPOL_INIT;
244 if (theta > 0) {
245 ptop = met->p[met->np - 1];
246 pbot = met->p[0];
247 do {
248 p0 = 0.5 * (ptop + pbot);
249 intpol_met_space_3d(met, met->t, p0, lons[ix], lats[iy],
250 &t0, ci, cw, 1);
251 if (THETA(p0, t0) > theta)
252 ptop = p0;
253 else
254 pbot = p0;
255 } while (fabs(ptop - pbot) > 1e-5);
256 }
257
258 /* Interpolate meteo data... */
259 INTPOL_SPACE_ALL(p0, lons[ix], lats[iy]);
260
261 /* Averaging... */
262 timem[iy * nx + ix] += met->time;
263 zm[iy * nx + ix] += z;
264 pm[iy * nx + ix] += p0;
265 tm[iy * nx + ix] += t;
266 um[iy * nx + ix] += u;
267 vm[iy * nx + ix] += v;
268 wm[iy * nx + ix] += w;
269 pvm[iy * nx + ix] += pv;
270 h2om[iy * nx + ix] += h2o;
271 o3m[iy * nx + ix] += o3;
272 lwcm[iy * nx + ix] += lwc;
273 rwcm[iy * nx + ix] += rwc;
274 iwcm[iy * nx + ix] += iwc;
275 swcm[iy * nx + ix] += swc;
276 ccm[iy * nx + ix] += cc;
277 psm[iy * nx + ix] += ps;
278 tsm[iy * nx + ix] += ts;
279 zsm[iy * nx + ix] += zs;
280 usm[iy * nx + ix] += us;
281 vsm[iy * nx + ix] += vs;
282 essm[iy * nx + ix] += ess;
283 nssm[iy * nx + ix] += nss;
284 shfm[iy * nx + ix] += shf;
285 lsmm[iy * nx + ix] += lsm;
286 sstm[iy * nx + ix] += sst;
287 pblm[iy * nx + ix] += pbl;
288 pctm[iy * nx + ix] += pct;
289 pcbm[iy * nx + ix] += pcb;
290 clm[iy * nx + ix] += cl;
291 if (isfinite(plfc) && isfinite(pel) && cape >= ctl.conv_cape
292 && (ctl.conv_cin <= 0 || cin < ctl.conv_cin)) {
293 plclm[iy * nx + ix] += plcl;
294 plfcm[iy * nx + ix] += plfc;
295 pelm[iy * nx + ix] += pel;
296 capem[iy * nx + ix] += cape;
297 cinm[iy * nx + ix] += cin;
298 npc[iy * nx + ix]++;
299 }
300 if (isfinite(pt)) {
301 ptm[iy * nx + ix] += pt;
302 ztm[iy * nx + ix] += zt;
303 ttm[iy * nx + ix] += tt;
304 h2otm[iy * nx + ix] += h2ot;
305 npt[iy * nx + ix]++;
306 }
307 o3cm[iy * nx + ix] += o3c;
308 hno3m[iy * nx + ix] += clim_zm(&clim->hno3, met->time, lats[iy], p0);
309 tnatm[iy * nx + ix] +=
310 nat_temperature(p0, h2o,
311 clim_zm(&clim->hno3, met->time, lats[iy], p0));
312 ohm[iy * nx + ix] +=
313 clim_oh(&ctl, clim, met->time, lons[ix], lats[iy], p0);
314 h2o2m[iy * nx + ix] += clim_zm(&clim->h2o2, met->time, lats[iy], p0);
315 ho2m[iy * nx + ix] += clim_zm(&clim->ho2, met->time, lats[iy], p0);
316 o1dm[iy * nx + ix] += clim_zm(&clim->o1d, met->time, lats[iy], p0);
317 rhm[iy * nx + ix] += RH(p0, t, h2o);
318 rhicem[iy * nx + ix] += RHICE(p0, t, h2o);
319 tdewm[iy * nx + ix] += TDEW(p0, h2o);
320 ticem[iy * nx + ix] += TICE(p0, h2o);
321 np[iy * nx + ix]++;
322 }
323 }
324
325 /* Create output file... */
326 LOG(1, "Write meteorological data file: %s", argv[2]);
327 if (!(out = fopen(argv[2], "w")))
328 ERRMSG("Cannot create file!");
329
330 /* Write header... */
331 MET_HEADER;
332
333 /* Write data... */
334 for (int iy = 0; iy < ny; iy++) {
335 fprintf(out, "\n");
336 for (int ix = 0; ix < nx; ix++)
337 fprintf(out,
338 "%.2f %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g"
339 " %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g"
340 " %g %g %g %g %g %g %g %g %g %g %g %g %g %g %d %d %d\n",
341 timem[iy * nx + ix] / np[iy * nx + ix],
342 Z(pm[iy * nx + ix] / np[iy * nx + ix]), lons[ix], lats[iy],
343 pm[iy * nx + ix] / np[iy * nx + ix],
344 tm[iy * nx + ix] / np[iy * nx + ix],
345 um[iy * nx + ix] / np[iy * nx + ix],
346 vm[iy * nx + ix] / np[iy * nx + ix],
347 wm[iy * nx + ix] / np[iy * nx + ix],
348 h2om[iy * nx + ix] / np[iy * nx + ix],
349 o3m[iy * nx + ix] / np[iy * nx + ix],
350 zm[iy * nx + ix] / np[iy * nx + ix],
351 pvm[iy * nx + ix] / np[iy * nx + ix],
352 psm[iy * nx + ix] / np[iy * nx + ix],
353 tsm[iy * nx + ix] / np[iy * nx + ix],
354 zsm[iy * nx + ix] / np[iy * nx + ix],
355 usm[iy * nx + ix] / np[iy * nx + ix],
356 vsm[iy * nx + ix] / np[iy * nx + ix],
357 essm[iy * nx + ix] / np[iy * nx + ix],
358 nssm[iy * nx + ix] / np[iy * nx + ix],
359 shfm[iy * nx + ix] / np[iy * nx + ix],
360 lsmm[iy * nx + ix] / np[iy * nx + ix],
361 sstm[iy * nx + ix] / np[iy * nx + ix],
362 ptm[iy * nx + ix] / npt[iy * nx + ix],
363 ztm[iy * nx + ix] / npt[iy * nx + ix],
364 ttm[iy * nx + ix] / npt[iy * nx + ix],
365 h2otm[iy * nx + ix] / npt[iy * nx + ix],
366 lwcm[iy * nx + ix] / np[iy * nx + ix],
367 rwcm[iy * nx + ix] / np[iy * nx + ix],
368 iwcm[iy * nx + ix] / np[iy * nx + ix],
369 swcm[iy * nx + ix] / np[iy * nx + ix],
370 ccm[iy * nx + ix] / np[iy * nx + ix],
371 clm[iy * nx + ix] / np[iy * nx + ix],
372 pctm[iy * nx + ix] / np[iy * nx + ix],
373 pcbm[iy * nx + ix] / np[iy * nx + ix],
374 plclm[iy * nx + ix] / npc[iy * nx + ix],
375 plfcm[iy * nx + ix] / npc[iy * nx + ix],
376 pelm[iy * nx + ix] / npc[iy * nx + ix],
377 capem[iy * nx + ix] / npc[iy * nx + ix],
378 cinm[iy * nx + ix] / npc[iy * nx + ix],
379 rhm[iy * nx + ix] / np[iy * nx + ix],
380 rhicem[iy * nx + ix] / np[iy * nx + ix],
381 tdewm[iy * nx + ix] / np[iy * nx + ix],
382 ticem[iy * nx + ix] / np[iy * nx + ix],
383 tnatm[iy * nx + ix] / np[iy * nx + ix],
384 hno3m[iy * nx + ix] / np[iy * nx + ix],
385 ohm[iy * nx + ix] / np[iy * nx + ix],
386 h2o2m[iy * nx + ix] / np[iy * nx + ix],
387 ho2m[iy * nx + ix] / np[iy * nx + ix],
388 o1dm[iy * nx + ix] / np[iy * nx + ix],
389 pblm[iy * nx + ix] / np[iy * nx + ix],
390 o3cm[iy * nx + ix] / np[iy * nx + ix], np[iy * nx + ix],
391 npt[iy * nx + ix], npc[iy * nx + ix]);
392 }
393
394 /* Close file... */
395 fclose(out);
396
397 /* Free... */
398 free(clim);
399 free(met);
400 free(dd);
401 free(timem);
402 free(psm);
403 free(tsm);
404 free(zsm);
405 free(usm);
406 free(vsm);
407 free(essm);
408 free(nssm);
409 free(shfm);
410 free(lsmm);
411 free(sstm);
412 free(pblm);
413 free(ptm);
414 free(pm);
415 free(tm);
416 free(um);
417 free(vm);
418 free(wm);
419 free(h2om);
420 free(h2otm);
421 free(o3m);
422 free(hno3m);
423 free(ohm);
424 free(h2o2m);
425 free(ho2m);
426 free(o1dm);
427 free(tdewm);
428 free(ticem);
429 free(tnatm);
430 free(lwcm);
431 free(rwcm);
432 free(iwcm);
433 free(swcm);
434 free(ccm);
435 free(zm);
436 free(pvm);
437 free(ztm);
438 free(ttm);
439 free(pctm);
440 free(pcbm);
441 free(clm);
442 free(plclm);
443 free(plfcm);
444 free(pelm);
445 free(capem);
446 free(cinm);
447 free(o3cm);
448 free(rhm);
449 free(rhicem);
450 free(np);
451 free(npc);
452 free(npt);
453
454 return EXIT_SUCCESS;
455}
NX
#define NX
Maximum number of longitudes.
Definition: met_map.c:32
NY
#define NY
Maximum number of latitudes.
Definition: met_map.c:35
clim_zm
double clim_zm(const clim_zm_t *zm, const double t, const double lat, const double p)
Interpolates monthly mean zonal mean climatological variables.
Definition: mptrac.c:405
nat_temperature
double nat_temperature(const double p, const double h2o, const double hno3)
Calculates the nitric acid trihydrate (NAT) temperature.
Definition: mptrac.c:7035
intpol_met_space_3d
void intpol_met_space_3d(const met_t *met, float array[EX][EY][EP], const double p, const double lon, const double lat, double *var, int *ci, double *cw, const int init)
Interpolates meteorological variables in 3D space.
Definition: mptrac.c:2301
scan_ctl
double scan_ctl(const char *filename, int argc, char *argv[], const char *varname, const int arridx, const char *defvalue, char *value)
Scans a control file or command-line arguments for a specified variable.
Definition: mptrac.c:11002
mptrac_read_clim
void mptrac_read_clim(const ctl_t *ctl, clim_t *clim)
Reads various climatological data and populates the given climatology structure.
Definition: mptrac.c:5501
clim_oh
double clim_oh(const ctl_t *ctl, const clim_t *clim, const double t, const double lon, const double lat, const double p)
Calculates the hydroxyl radical (OH) concentration from climatology data, with an optional diurnal co...
Definition: mptrac.c:89
mptrac_read_met
int mptrac_read_met(const char *filename, const ctl_t *ctl, const clim_t *clim, met_t *met, dd_t *dd)
Reads meteorological data from a file, supporting multiple formats and MPI broadcasting.
Definition: mptrac.c:6467
mptrac_read_ctl
void mptrac_read_ctl(const char *filename, int argc, char *argv[], ctl_t *ctl)
Reads control parameters from a configuration file and populates the given structure.
Definition: mptrac.c:5561
INTPOL_SPACE_ALL
#define INTPOL_SPACE_ALL(p, lon, lat)
Interpolate multiple meteorological variables in space.
Definition: mptrac.h:933
INTPOL_INIT
#define INTPOL_INIT
Initialize arrays for interpolation.
Definition: mptrac.h:883
ERRMSG
#define ERRMSG(...)
Print an error message with contextual information and terminate the program.
Definition: mptrac.h:2115
USAGE
#define USAGE
Print usage information on -h or --help.
Definition: mptrac.h:1922
Z
#define Z(p)
Convert pressure to altitude.
Definition: mptrac.h:1952
P
#define P(z)
Compute pressure at given altitude.
Definition: mptrac.h:1493
THETA
#define THETA(p, t)
Compute potential temperature.
Definition: mptrac.h:1833
MET_HEADER
#define MET_HEADER
Write header for meteorological data file.
Definition: mptrac.h:1101
TICE
#define TICE(p, h2o)
Calculate frost point temperature (WMO, 2018).
Definition: mptrac.h:1809
RHICE
#define RHICE(p, t, h2o)
Compute relative humidity over ice.
Definition: mptrac.h:1645
ALLOC
#define ALLOC(ptr, type, n)
Allocate memory for a pointer with error handling.
Definition: mptrac.h:466
RH
#define RH(p, t, h2o)
Compute relative humidity over water.
Definition: mptrac.h:1615
LOG
#define LOG(level,...)
Print a log message with a specified logging level.
Definition: mptrac.h:2045
TDEW
#define TDEW(p, h2o)
Calculate dew point temperature.
Definition: mptrac.h:1784
clim_t
Climatological data.
Definition: mptrac.h:3452
clim_t::ho2
clim_zm_t ho2
HO2 zonal means.
Definition: mptrac.h:3482
clim_t::hno3
clim_zm_t hno3
HNO3 zonal means.
Definition: mptrac.h:3473
clim_t::o1d
clim_zm_t o1d
O(1D) zonal means.
Definition: mptrac.h:3485
clim_t::h2o2
clim_zm_t h2o2
H2O2 zonal means.
Definition: mptrac.h:3479
ctl_t
Control parameters.
Definition: mptrac.h:2203
ctl_t::conv_cape
double conv_cape
CAPE threshold for convection module [J/kg].
Definition: mptrac.h:2778
ctl_t::conv_cin
double conv_cin
CIN threshold for convection module [J/kg].
Definition: mptrac.h:2781
dd_t
Domain decomposition data structure.
Definition: mptrac.h:3685
met_t
Meteo data structure.
Definition: mptrac.h:3511
met_t::nx
int nx
Number of longitudes.
Definition: mptrac.h:3517
met_t::ny
int ny
Number of latitudes.
Definition: mptrac.h:3520
met_t::np
int np
Number of pressure levels.
Definition: mptrac.h:3523
met_t::t
float t[EX][EY][EP]
Temperature [K].
Definition: mptrac.h:3625
met_t::lon
double lon[EX]
Longitudes [deg].
Definition: mptrac.h:3529
met_t::time
double time
Time [s].
Definition: mptrac.h:3514
met_t::lat
double lat[EY]
Latitudes [deg].
Definition: mptrac.h:3532
met_t::p
double p[EP]
Pressure levels [hPa].
Definition: mptrac.h:3535
Here is the call graph for this function:
Generated by doxygen 1.9.4