1 /*
2 * (C) Copyright 2005- ECMWF.
3 *
4 * This software is licensed under the terms of the Apache Licence Version 2.0
5 * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
6 *
7 * In applying this licence, ECMWF does not waive the privileges and immunities granted to it by
8 * virtue of its status as an intergovernmental organisation nor does it submit to any jurisdiction.
9 */
10
11 #include "grib_api_internal.h"
12 #include <math.h>
13
14 /*
15 This is used by make_class.pl
16
17 START_CLASS_DEF
18 CLASS = iterator
19 SUPER = grib_iterator_class_gen
20 IMPLEMENTS = previous;next
21 IMPLEMENTS = init;destroy
22 MEMBERS = double *las
23 MEMBERS = double *los
24 MEMBERS = long Ni
25 MEMBERS = long Nj
26 MEMBERS = long iScansNegatively
27 MEMBERS = long isRotated
28 MEMBERS = double angleOfRotation
29 MEMBERS = double southPoleLat
30 MEMBERS = double southPoleLon
31 MEMBERS = long jPointsAreConsecutive
32 MEMBERS = long disableUnrotate
33 END_CLASS_DEF
34
35 */
36
37 /* START_CLASS_IMP */
38
39 /*
40
41 Don't edit anything between START_CLASS_IMP and END_CLASS_IMP
42 Instead edit values between START_CLASS_DEF and END_CLASS_DEF
43 or edit "iterator.class" and rerun ./make_class.pl
44
45 */
46
47
48 static void init_class(grib_iterator_class*);
49
50 static int init(grib_iterator* i, grib_handle*, grib_arguments*);
51 static int next(grib_iterator* i, double* lat, double* lon, double* val);
52 static int previous(grib_iterator* ei, double* lat, double* lon, double* val);
53 static int destroy(grib_iterator* i);
54
55
56 typedef struct grib_iterator_regular
57 {
58 grib_iterator it;
59 /* Members defined in gen */
60 long carg;
61 const char* missingValue;
62 /* Members defined in regular */
63 double* las;
64 double* los;
65 long Ni;
66 long Nj;
67 long iScansNegatively;
68 long isRotated;
69 double angleOfRotation;
70 double southPoleLat;
71 double southPoleLon;
72 long jPointsAreConsecutive;
73 long disableUnrotate;
74 } grib_iterator_regular;
75
76 extern grib_iterator_class* grib_iterator_class_gen;
77
78 static grib_iterator_class _grib_iterator_class_regular = {
79 &grib_iterator_class_gen, /* super */
80 "regular", /* name */
81 sizeof(grib_iterator_regular), /* size of instance */
82 0, /* inited */
83 &init_class, /* init_class */
84 &init, /* constructor */
85 &destroy, /* destructor */
86 &next, /* Next Value */
87 &previous, /* Previous Value */
88 0, /* Reset the counter */
89 0, /* has next values */
90 };
91
92 grib_iterator_class* grib_iterator_class_regular = &_grib_iterator_class_regular;
93
94
init_class(grib_iterator_class * c)95 static void init_class(grib_iterator_class* c)
96 {
97 c->reset = (*(c->super))->reset;
98 c->has_next = (*(c->super))->has_next;
99 }
100 /* END_CLASS_IMP */
101
102
next(grib_iterator * i,double * lat,double * lon,double * val)103 static int next(grib_iterator* i, double* lat, double* lon, double* val)
104 {
105 grib_iterator_regular* self = (grib_iterator_regular*)i;
106
107 if ((long)i->e >= (long)(i->nv - 1))
108 return 0;
109
110 i->e++;
111
112 *lat = self->las[(long)floor(i->e / self->Ni)];
113 *lon = self->los[(long)i->e % self->Ni];
114 *val = i->data[i->e];
115
116 return 1;
117 }
118
previous(grib_iterator * i,double * lat,double * lon,double * val)119 static int previous(grib_iterator* i, double* lat, double* lon, double* val)
120 {
121 grib_iterator_regular* self = (grib_iterator_regular*)i;
122
123 if (i->e < 0)
124 return 0;
125 *lat = self->las[(long)floor(i->e / self->Ni)];
126 *lon = self->los[i->e % self->Ni];
127 *val = i->data[i->e];
128 i->e--;
129
130 return 1;
131 }
132
destroy(grib_iterator * i)133 static int destroy(grib_iterator* i)
134 {
135 grib_iterator_regular* self = (grib_iterator_regular*)i;
136 const grib_context* c = i->h->context;
137 grib_context_free(c, self->las);
138 grib_context_free(c, self->los);
139 return GRIB_SUCCESS;
140 }
141
init(grib_iterator * i,grib_handle * h,grib_arguments * args)142 static int init(grib_iterator* i, grib_handle* h, grib_arguments* args)
143 {
144 grib_iterator_regular* self = (grib_iterator_regular*)i;
145 int ret = GRIB_SUCCESS;
146
147 long Ni; /* Number of points along a parallel = Nx */
148 long Nj; /* Number of points along a meridian = Ny */
149 double idir, lon1, lon2;
150 long loi;
151
152 const char* s_lon1 = grib_arguments_get_name(h, args, self->carg++);
153 const char* s_idir = grib_arguments_get_name(h, args, self->carg++);
154 const char* s_Ni = grib_arguments_get_name(h, args, self->carg++);
155 const char* s_Nj = grib_arguments_get_name(h, args, self->carg++);
156 const char* s_iScansNeg = grib_arguments_get_name(h, args, self->carg++);
157
158 if ((ret = grib_get_double_internal(h, s_lon1, &lon1)))
159 return ret;
160 if ((ret = grib_get_double_internal(h, "longitudeOfLastGridPointInDegrees", &lon2)))
161 return ret;
162 if ((ret = grib_get_double_internal(h, s_idir, &idir)))
163 return ret;
164 if ((ret = grib_get_long_internal(h, s_Ni, &Ni)))
165 return ret;
166 if (grib_is_missing(h, s_Ni, &ret) && ret == GRIB_SUCCESS) {
167 grib_context_log(h->context, GRIB_LOG_ERROR, "Key %s cannot be 'missing' for a regular grid!", s_Ni);
168 return GRIB_WRONG_GRID;
169 }
170
171 if ((ret = grib_get_long_internal(h, s_Nj, &Nj)))
172 return ret;
173 if (grib_is_missing(h, s_Nj, &ret) && ret == GRIB_SUCCESS) {
174 grib_context_log(h->context, GRIB_LOG_ERROR, "Key %s cannot be 'missing' for a regular grid!", s_Nj);
175 return GRIB_WRONG_GRID;
176 }
177
178 if ((ret = grib_get_long_internal(h, s_iScansNeg, &self->iScansNegatively)))
179 return ret;
180
181 /* GRIB-801: Careful of case with a single point! Ni==1 */
182 if (Ni > 1) {
183 /* Note: If first and last longitudes are equal I assume you wanna go round the globe */
184 if (self->iScansNegatively) {
185 if (lon1 > lon2) {
186 idir = (lon1 - lon2) / (Ni - 1);
187 }
188 else {
189 idir = (lon1 + 360.0 - lon2) / (Ni - 1);
190 }
191 }
192 else {
193 if (lon2 > lon1) {
194 idir = (lon2 - lon1) / (Ni - 1);
195 }
196 else {
197 idir = (lon2 + 360.0 - lon1) / (Ni - 1);
198 }
199 }
200 }
201 if (self->iScansNegatively) {
202 idir = -idir;
203 }
204 else {
205 if (lon1 + (Ni - 2) * idir > 360)
206 lon1 -= 360;
207 /*See ECC-704, GRIB-396*/
208 /*else if ( (lon1+(Ni-1)*idir)-360 > epsilon ){
209 idir=360.0/(float)Ni;
210 }*/
211 }
212
213 self->Ni = Ni;
214 self->Nj = Nj;
215
216 self->las = (double*)grib_context_malloc(h->context, Nj * sizeof(double));
217 self->los = (double*)grib_context_malloc(h->context, Ni * sizeof(double));
218
219 for (loi = 0; loi < Ni; loi++) {
220 self->los[loi] = lon1;
221 lon1 += idir;
222 }
223
224 return ret;
225 }
226