MPTRAC
Functions
atm_select.c File Reference

Extract subsets of air parcels from atmospheric data files. More...

#include "mptrac.h"

Go to the source code of this file.

Functions

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

Detailed Description

Extract subsets of air parcels from atmospheric data files.

Definition in file atm_select.c.

Function Documentation

◆ main()

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

Definition at line 27 of file atm_select.c.

29 {
30
31 ctl_t ctl;
32
33 atm_t *atm, *atm2;
34
35 /* Allocate... */
36 ALLOC(atm, atm_t, 1);
37 ALLOC(atm2, atm_t, 1);
38
39 /* Check arguments... */
40 if (argc < 4)
41 ERRMSG("Give parameters: <ctl> <atm_select> <atm1> [<atm2> ...]");
42
43 /* Read control parameters... */
44 read_ctl(argv[1], argc, argv, &ctl);
45 const int stride =
46 (int) scan_ctl(argv[1], argc, argv, "SELECT_STRIDE", -1, "1", NULL);
47 const int idx0 =
48 (int) scan_ctl(argv[1], argc, argv, "SELECT_IDX0", -1, "-999", NULL);
49 const int idx1 =
50 (int) scan_ctl(argv[1], argc, argv, "SELECT_IDX1", -1, "-999", NULL);
51 int ip0 =
52 (int) scan_ctl(argv[1], argc, argv, "SELECT_IP0", -1, "-999", NULL);
53 int ip1 =
54 (int) scan_ctl(argv[1], argc, argv, "SELECT_IP1", -1, "-999", NULL);
55 const double t0 = scan_ctl(argv[1], argc, argv, "SELECT_T0", -1, "0", NULL);
56 const double t1 = scan_ctl(argv[1], argc, argv, "SELECT_T1", -1, "0", NULL);
57 const double p0 =
58 P(scan_ctl(argv[1], argc, argv, "SELECT_Z0", -1, "0", NULL));
59 const double p1 =
60 P(scan_ctl(argv[1], argc, argv, "SELECT_Z1", -1, "0", NULL));
61 const double theta0 =
62 scan_ctl(argv[1], argc, argv, "SELECT_THETA0", -1, "0", NULL);
63 const double theta1 =
64 scan_ctl(argv[1], argc, argv, "SELECT_THETA1", -1, "0", NULL);
65 const double lon0 =
66 scan_ctl(argv[1], argc, argv, "SELECT_LON0", -1, "0", NULL);
67 const double lon1 =
68 scan_ctl(argv[1], argc, argv, "SELECT_LON1", -1, "0", NULL);
69 const double lat0 =
70 scan_ctl(argv[1], argc, argv, "SELECT_LAT0", -1, "0", NULL);
71 const double lat1 =
72 scan_ctl(argv[1], argc, argv, "SELECT_LAT1", -1, "0", NULL);
73 const double r0 = scan_ctl(argv[1], argc, argv, "SELECT_R0", -1, "0", NULL);
74 const double r1 = scan_ctl(argv[1], argc, argv, "SELECT_R1", -1, "0", NULL);
75 const double rlon =
76 scan_ctl(argv[1], argc, argv, "SELECT_RLON", -1, "0", NULL);
77 const double rlat =
78 scan_ctl(argv[1], argc, argv, "SELECT_RLAT", -1, "0", NULL);
79
80 /* Get Cartesian coordinates... */
81 double x0[3], x1[3];
82 geo2cart(0, rlon, rlat, x0);
83
84 /* Loop over files... */
85 for (int f = 3; f < argc; f++) {
86
87 /* Read atmopheric data... */
88 if (!read_atm(argv[f], &ctl, atm))
89 continue;
90
91 /* Adjust range of air parcels... */
92 if (ip0 < 0)
93 ip0 = 0;
94 ip0 = MIN(ip0, atm->np - 1);
95 if (ip1 < 0)
96 ip1 = atm->np - 1;
97 ip1 = MIN(ip1, atm->np - 1);
98 if (ip1 < ip0)
99 ip1 = ip0;
100
101 /* Loop over air parcels... */
102 for (int ip = ip0; ip <= ip1; ip += stride) {
103
104 /* Check air parcel index... */
105 if (ctl.qnt_idx >= 0 && idx0 >= 0 && idx1 >= 0)
106 if (atm->q[ctl.qnt_idx][ip] < idx0 || atm->q[ctl.qnt_idx][ip] > idx1)
107 continue;
108
109 /* Check time... */
110 if (t0 != t1)
111 if ((t1 > t0 && (atm->time[ip] < t0 || atm->time[ip] > t1))
112 || (t1 < t0 && (atm->time[ip] < t0 && atm->time[ip] > t1)))
113 continue;
114
115 /* Check vertical distance... */
116 if (p0 != p1)
117 if ((p0 > p1 && (atm->p[ip] > p0 || atm->p[ip] < p1))
118 || (p0 < p1 && (atm->p[ip] > p0 && atm->p[ip] < p1)))
119 continue;
120
121 /* Check potential temperature... */
122 if (theta0 != theta1) {
123 double theta;
124 if (ctl.qnt_theta >= 0)
125 theta = atm->q[ctl.qnt_theta][ip];
126 else if (ctl.qnt_t >= 0)
127 theta = THETA(atm->p[ip], atm->q[ctl.qnt_t][ip]);
128 else
129 ERRMSG
130 ("Filtering requires temperature or potential temperature data!");
131 if ((theta1 > theta0 && (theta < theta0 || theta > theta1))
132 || (theta1 < theta0 && (theta < theta0 && theta > theta1)))
133 continue;
134 }
135
136 /* Check longitude... */
137 if (lon0 != lon1)
138 if ((lon1 > lon0 && (atm->lon[ip] < lon0 || atm->lon[ip] > lon1))
139 || (lon1 < lon0 && (atm->lon[ip] < lon0 && atm->lon[ip] > lon1)))
140 continue;
141
142 /* Check latitude... */
143 if (lat0 != lat1)
144 if ((lat1 > lat0 && (atm->lat[ip] < lat0 || atm->lat[ip] > lat1))
145 || (lat1 < lat0 && (atm->lat[ip] < lat0 && atm->lat[ip] > lat1)))
146 continue;
147
148 /* Check horizontal distace... */
149 if (r0 != r1) {
150 geo2cart(0, atm->lon[ip], atm->lat[ip], x1);
151 double r = DIST(x0, x1);
152 if ((r1 > r0 && (r < r0 || r > r1))
153 || (r1 < r0 && (r < r0 && r > r1)))
154 continue;
155 }
156
157 /* Copy data... */
158 atm2->time[atm2->np] = atm->time[ip];
159 atm2->p[atm2->np] = atm->p[ip];
160 atm2->lon[atm2->np] = atm->lon[ip];
161 atm2->lat[atm2->np] = atm->lat[ip];
162 for (int iq = 0; iq < ctl.nq; iq++)
163 atm2->q[iq][atm2->np] = atm->q[iq][ip];
164 if ((++atm2->np) > NP)
165 ERRMSG("Too many air parcels!");
166 }
167 }
168
169 /* Close file... */
170 write_atm(argv[2], &ctl, atm2, 0);
171
172 /* Free... */
173 free(atm);
174 free(atm2);
175
176 return EXIT_SUCCESS;
177}
write_atm
void write_atm(const char *filename, const ctl_t *ctl, const atm_t *atm, const double t)
Writes air parcel data to a file in various formats.
Definition: mptrac.c:8874
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:8503
read_ctl
void 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:5156
read_atm
int read_atm(const char *filename, const ctl_t *ctl, atm_t *atm)
Reads air parcel data from a specified file into the given atmospheric structure.
Definition: mptrac.c:4566
geo2cart
void geo2cart(const double z, const double lon, const double lat, double *x)
Converts geographic coordinates (longitude, latitude, altitude) to Cartesian coordinates.
Definition: mptrac.c:985
MIN
#define MIN(a, b)
Macro to determine the minimum of two values.
Definition: mptrac.h:987
ERRMSG
#define ERRMSG(...)
Print an error message with contextual information and terminate the program.
Definition: mptrac.h:1916
P
#define P(z)
Compute pressure at given altitude.
Definition: mptrac.h:1304
THETA
#define THETA(p, t)
Compute potential temperature.
Definition: mptrac.h:1644
ALLOC
#define ALLOC(ptr, type, n)
Allocate memory for a pointer with error handling.
Definition: mptrac.h:349
NP
#define NP
Maximum number of atmospheric data points.
Definition: mptrac.h:246
DIST
#define DIST(a, b)
Calculate the distance between two points in Cartesian coordinates.
Definition: mptrac.h:578
atm_t
Air parcel data.
Definition: mptrac.h:3147
atm_t::time
double time[NP]
Time [s].
Definition: mptrac.h:3153
atm_t::lat
double lat[NP]
Latitude [deg].
Definition: mptrac.h:3162
atm_t::lon
double lon[NP]
Longitude [deg].
Definition: mptrac.h:3159
atm_t::np
int np
Number of air parcels.
Definition: mptrac.h:3150
atm_t::q
double q[NQ][NP]
Quantity data (for various, user-defined attributes).
Definition: mptrac.h:3165
atm_t::p
double p[NP]
Pressure [hPa].
Definition: mptrac.h:3156
ctl_t
Control parameters.
Definition: mptrac.h:2170
ctl_t::qnt_theta
int qnt_theta
Quantity array index for potential temperature.
Definition: mptrac.h:2381
ctl_t::qnt_t
int qnt_t
Quantity array index for temperature.
Definition: mptrac.h:2264
ctl_t::qnt_idx
int qnt_idx
Quantity array index for air parcel IDs.
Definition: mptrac.h:2192
ctl_t::nq
int nq
Number of quantities.
Definition: mptrac.h:2177
Here is the call graph for this function:
Generated on Wed Jan 22 2025 18:02:19 for MPTRAC by doxygen 1.9.4