1 /*
2   This file is part of CDO. CDO is a collection of Operators to manipulate and analyse Climate model Data.
3 
4   Copyright (C) 2006 Brockmann Consult
5 
6   Author: Ralf Quast
7 
8 */
9 #ifndef ECAUTIL_H_
10 #define ECAUTIL_H_
11 
12 #include "field.h"
13 
14 /**
15  * Computes the day-of-year correspnding a given Gregorian date.
16  *
17  * @param date a Gregorian date in the form YYYYMMDD
18  *
19  * @return the day-of-year
20  */
21 unsigned long day_of_year(int date);
22 
23 /**
24  * Counts the number of nonmissing values in a field. The result
25  * of the operation is computed according to the following rules:
26  *
27  * field1  field2  result
28  * a       b       a + 1
29  * a       miss    a
30  * miss    b       1
31  * miss    miss    0
32  *
33  * @param field1 the 1st input field, also holds the result
34  * @param field2 the 2nd input field
35  */
36 void vfarnum(Field &field1, const Field &field2);
37 
38 /**
39  * Counts the number of consecutive nonmissing values in a field.
40  * The result of the operation is computed according to the following
41  * rules:
42  *
43  * field1  field2  result
44  * a       b       a + 1
45  * a       miss    0
46  * miss    b       1
47  * miss    miss    0
48  *
49  * @param field1 the 1st input field, also holds the result
50  * @param field2 the 2nd input field
51  */
52 void vfarnum2(Field &field1, const Field &field2);
53 
54 /**
55  * Counts the number of values in series of at least n consecutive
56  * nonmissing values. The result of the operation is computed according
57  * to the following rules:
58  *
59  * field1  field2  result
60  * a       b       b < n ? a : b > n ? a + 1 : a + n
61  * a       miss    a
62  * miss    b       b < n ? 0 : b
63  * miss    miss    0
64  *
65  * @param field1 the 1st input field, also holds the result
66  * @param field2 the 2nd input field
67  * @param n      the number of consecutive values, must be an exact
68  *               mathematical integer
69  */
70 void vfarnum3(Field &field1, const Field &field2, double n);
71 
72 /**
73  * Selects field elements according to a given mask. The result of
74  * the operation is computed according to the following rules:
75  *
76  * field1  field2  result
77  * a       b       b != 0.0 ? a : miss
78  * a       miss    miss
79  * miss    b       miss
80  * miss    miss    miss
81  *
82  * @param field1  the input field, also holds the result
83  * @param field2  the mask
84  */
85 void vfarsel(Field &field1, const Field &field2);
86 
87 /**
88  * Selects all field elements that are less than or equal to the
89  * corresponding element of a reference field. The result of the
90  * operation is computed according to the following rules:
91  *
92  * field1  field2  result
93  * a       b       a <= b ? a : miss
94  * a       miss    miss
95  * miss    b       miss
96  * miss    miss    miss
97  *
98  * @param field1 the input field, also holds the result
99  * @param field2 the reference field
100  */
101 void vfarselle(Field &field1, const Field &field2);
102 
103 /**
104  * Selects all field elements that are less than the
105  * corresponding element of a reference field. The result of the
106  * operation is computed according to the following rules:
107  *
108  * field1  field2  result
109  * a       b       a < b ? a : miss
110  * a       miss    miss
111  * miss    b       miss
112  * miss    miss    miss
113  *
114  * @param field1 the input field, also holds the result
115  * @param field2 the reference field
116  */
117 void vfarsellt(Field &field1, const Field &field2);
118 
119 /**
120  * Selects all field elements that are greater than or equal to
121  * the corresponding element of a reference field. The result of
122  * the operation is computed according to the following rules:
123  *
124  * field1  field2  result
125  * a       b       a >= b ? a : miss
126  * a       miss    miss
127  * miss    b       miss
128  * miss    miss    miss
129  *
130  * @param field1 the input field, also holds the result
131  * @param field2 the reference field
132  */
133 void vfarselge(Field &field1, const Field &field2);
134 
135 /**
136  * Selects all field elements that are greater than the
137  * corresponding element of a reference field. The result of the
138  * operation is computed according to the following rules:
139  *
140  * field1  field2  result
141  * a       b       a > b ? a : miss
142  * a       miss    miss
143  * miss    b       miss
144  * miss    miss    miss
145  *
146  * @param field1 the input field, also holds the result
147  * @param field2 the reference field
148  */
149 void vfarselgt(Field &field1, const Field &field2);
150 
151 /**
152  * Selects all field elements that are equal to the
153  * corresponding element of a reference field. The result of the
154  * operation is computed according to the following rules:
155  *
156  * field1  field2  result
157  * a       b       a == b ? a : miss
158  * a       miss    miss
159  * miss    b       miss
160  * miss    miss    miss
161  *
162  * @param field1 the input field, also holds the result
163  * @param field2 the reference field
164  */
165 void vfarseleq(Field &field1, const Field &field2);
166 
167 /**
168  * Selects all field elements that are not equal to the
169  * corresponding element of a reference field. The result of the
170  * operation is computed according to the following rules:
171  *
172  * field1  field2  result
173  * a       b       a != b ? a : miss
174  * a       miss    miss
175  * miss    b       miss
176  * miss    miss    miss
177  *
178  * @param field1 the input field, also holds the result
179  * @param field2 the reference field
180  */
181 void vfarselne(Field &field1, const Field &field2);
182 
183 /**
184  * Selects all field elements that are less than or equal to a
185  * certain reference value. The result of the operation is computed
186  * according to the following rules:
187  *
188  * field  c      result
189  * a      c      a <= c ? a : miss
190  * a      miss   miss
191  * miss   c      miss
192  * miss   miss   miss
193  *
194  * @param field the input field, also holds the result
195  * @param c     the reference value
196  */
197 void vfarsellec(Field &field, double c);
198 
199 /**
200  * Selects all field elements that are less a
201  * certain reference value. The result of the operation is computed
202  * according to the following rules:
203  *
204  * field  c      result
205  * a      c      a < c ? a : miss
206  * a      miss   miss
207  * miss   c      miss
208  * miss   miss   miss
209  *
210  * @param field the input field, also holds the result
211  * @param c     the reference value
212  */
213 void vfarselltc(Field &field, double c);
214 
215 /**
216  * Selects all field elements that are greater than or equal to a
217  * certain reference value. The result of the operation is computed
218  * according to the following rules:
219  *
220  * field  c      result
221  * a      c      a >= c ? a : miss
222  * a      miss   miss
223  * miss   c      miss
224  * miss   miss   miss
225  *
226  * @param field the input field, also holds the result
227  * @param c     the reference value
228  */
229 void vfarselgec(Field &field, double c);
230 
231 /**
232  * Selects all field elements that are greater than a
233  * certain reference value. The result of the operation is computed
234  * according to the following rules:
235  *
236  * field  c      result
237  * a      c      a > c ? a : miss
238  * a      miss   miss
239  * miss   c      miss
240  * miss   miss   miss
241  *
242  * @param field the input field, also holds the result
243  * @param c     the reference value
244  */
245 void vfarselgtc(Field &field, double c);
246 
247 /**
248  * Selects all field elements that are equal to a
249  * certain reference value. The result of the operation is computed
250  * according to the following rules:
251  *
252  * field  c      result
253  * a      c      a == c ? a : miss
254  * a      miss   miss
255  * miss   c      miss
256  * miss   miss   miss
257  *
258  * @param field the input field, also holds the result
259  * @param c     the reference value
260  */
261 void vfarseleqc(Field &field, double c);
262 
263 /**
264  * Selects all field elements that are not equal to a
265  * certain reference value. The result of the operation is computed
266  * according to the following rules:
267  *
268  * field  c      result
269  * a      c      a != c ? a : miss
270  * a      miss   miss
271  * miss   c      miss
272  * miss   miss   miss
273  *
274  * @param field the input field, also holds the result
275  * @param c     the reference value
276  */
277 void vfarselnec(Field &field, double c);
278 
279 /**
280  * reset the fields real values to the missval for all levels
281  *
282  * @param field     list of fields: 0 is index of the current values, 1 hold
283  *                  the values of the previous year
284  * @param nlevels   number of available levels
285  * @param gridsize  number of grid points
286  * @param yvals     list of latitudes
287  * @param onlyNorth boolean for processing only the norther hemisphere
288  */
289 void update_hist(FieldVector2D &field, int nlevels, size_t gridsize, const std::vector<double> &yvals, bool onlyNorth);
290 
291 /*
292  * Compute the Gsl and its starting day
293  *
294  * @param int nlevels
295  * @param size_t gridsize
296  * @param double *yvals = array of latitudes
297  * @param int ysize = number of gridpoints in lat-direction
298  * @param double missval
299  * @param int ovdate = the last year, which has been fully processed
300  * @param Field *startDate
301  * @param Field *endDate
302  * @param Field *startDateWithHist[2]
303  * @param Field *endDateWithHist[2]
304  * @param Field *gslDuration
305  * @param Field *gslFirstDay
306  * @param bool useCurrentYear = if true, only measurements of the current year
307  *                             (index 0) are used for computation, i.e. that
308  *                             gsl can only be computed for the northern
309  *                             hemisphere (see definition of GSL: EcaGsl()
310  */
311 void compute_gsl(int nlevels, size_t gridsize, std::vector<double> &yvals, double missval, FieldVector2D &startDateWithHist,
312                  FieldVector2D &endDateWithHist, FieldVector &gslDuration, FieldVector &gslFirstDay, bool useCurrentYear);
313 
314 /*
315  * Adjust the endDates found in the current year:
316  * if a starting date for gsl could be found, but no ending date, the end
317  * should be the last day of the corresponding year for norther and June, 30th
318  * for southern hemisphere
319  */
320 void adjust_end_date(int nlevels, size_t gridsize, const std::vector<double> &yvals, double missval, int64_t ovdate,
321                      const FieldVector2D &startDateWithHist, FieldVector2D &endDateWithHist);
322 /*
323  * Calculates the mid of the year or month
324  */
325 void define_mid_of_time(int frequency, int taxisID, int year, int month, int MaxMonths);
326 /*
327  * Write GSL related fields to an output stream
328  */
329 void write_gsl_stream(CdoStreamID ostreamID, int otaxisID, int otsID, int ovarID1, int ovarID2, int ivlistID1, int first_var_id,
330                       FieldVector &gslDuration, FieldVector &gslFirstDay, int64_t vdate, int vtime, int nlevels);
331 #endif /*ECAUTIL_H_*/
332