bes  Updated for version 3.20.13
HDFEOS2ArraySwathGeoMultiDimMapField.cc
1 // Retrieves the latitude and longitude of the HDF-EOS2 Swath without using dimension maps
3 //
4 // Authors: MuQun Yang <myang6@hdfgroup.org>
5 // Copyright (c) 2010-2020 The HDF Group
7 //For the swath using the multiple dimension maps,
8 // Latitude/longitude will be interpolated accordingly.
9 
10 #ifdef USE_HDFEOS2_LIB
11 
12 #include "HDFEOS2ArraySwathGeoMultiDimMapField.h"
13 #include <iostream>
14 #include <sstream>
15 #include <cassert>
16 #include <libdap/debug.h>
17 #include <libdap/InternalErr.h>
18 #include "BESDebug.h"
19 #include "HDF4RequestHandler.h"
20 
21 using namespace std;
22 using namespace libdap;
23 #define SIGNED_BYTE_TO_INT32 1
24 
25 bool
26 HDFEOS2ArraySwathGeoMultiDimMapField::read ()
27 {
28 
29  BESDEBUG("h4","Coming to HDFEOS2ArraySwathGeoMultiDimMapField read "<<endl);
30 
31  if(length() == 0)
32  return true;
33 
34  if(rank !=2)
35  throw InternalErr (__FILE__, __LINE__, "The field rank must be 2 for swath multi-dimension map reading.");
36 
37  bool check_pass_fileid_key = HDF4RequestHandler::get_pass_fileid();
38 
39  // Declare offset, count and step
40  vector<int>offset;
41  offset.resize(rank);
42  vector<int>count;
43  count.resize(rank);
44  vector<int>step;
45  step.resize(rank);
46 
47  // Obtain offset,step and count from the client expression constraint
48  int nelms = format_constraint (offset.data(), step.data(), count.data());
49 
50  // Just declare offset,count and step in the int32 type.
51  vector<int32>offset32;
52  offset32.resize(rank);
53  vector<int32>count32;
54  count32.resize(rank);
55  vector<int32>step32;
56  step32.resize(rank);
57 
58  // Just obtain the offset,count and step in the datatype of int32.
59  for (int i = 0; i < rank; i++) {
60  offset32[i] = (int32) offset[i];
61  count32[i] = (int32) count[i];
62  step32[i] = (int32) step[i];
63  }
64 
65  int32 (*openfunc) (char *, intn);
66  int32 (*attachfunc) (int32, char *);
67  intn (*detachfunc) (int32);
68  intn (*fieldinfofunc) (int32, char *, int32 *, int32 *, int32 *, char *);
69  intn (*readfieldfunc) (int32, char *, int32 *, int32 *, int32 *, void *);
70 
71  // Define function pointers to handle the swath
72  string datasetname;
73  openfunc = SWopen;
74  attachfunc = SWattach;
75  detachfunc = SWdetach;
76  fieldinfofunc = SWfieldinfo;
77  readfieldfunc = SWreadfield;
78  datasetname = swathname;
79 
80  // Declare dimension map size, offset and inc
81  vector<int>dm_dimsizes;
82  dm_dimsizes.resize(rank);
83  vector<int>dm_offsets;
84  dm_offsets.resize(rank);
85  vector<int>dm_incs;
86  dm_incs.resize(rank);
87  bool no_interpolation = true;
88 
89  if(dim0inc != 0 || dim1inc !=0 || dim0offset != 0 || dim1offset != 0) {
90  dm_dimsizes[0] = dim0size;
91  dm_dimsizes[1] = dim1size;
92  dm_offsets[0] = dim0offset;
93  dm_offsets[1] = dim1offset;
94  dm_incs[0] = dim0inc;
95  dm_incs[1] = dim1inc;
96  no_interpolation = false;
97  }
98 
99  // We may eventually combine the following code with other code, so
100  // we don't add many comments from here to the end of the file.
101  // The jira ticket about combining code is HFRHANDLER-166.
102  int32 sfid = -1;
103  int32 swathid = -1;
104 
105  if(false == check_pass_fileid_key) {
106  sfid = openfunc (const_cast < char *>(filename.c_str ()), DFACC_READ);
107  if (sfid < 0) {
108  ostringstream eherr;
109  eherr << "File " << filename.c_str () << " cannot be open.";
110  throw InternalErr (__FILE__, __LINE__, eherr.str ());
111  }
112  }
113  else
114  sfid = swathfd;
115 
116 #if 0
117 cerr<<"swath name is "<<datasetname <<endl;
118 cerr<<"swath rank is "<<rank <<endl;
119 cerr<<"swath field name is "<<fieldname <<endl;
120 cerr<<"swath field new name is "<<name() <<endl;
121 
122 cerr<<"datadim0size "<<dim0size <<endl;
123 cerr<<"datadim1size "<<dim1size <<endl;
124 cerr<<"datadim0offset "<<dim0offset <<endl;
125 cerr<<"datadim1offset "<<dim1offset <<endl;
126 cerr<<"datadim0inc "<<dim0inc <<endl;
127 cerr<<"datadim1inc "<<dim1inc <<endl;
128 #endif
129 
130  swathid = attachfunc (sfid, const_cast < char *>(datasetname.c_str ()));
131 
132  if (swathid < 0) {
133  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
134  ostringstream eherr;
135  eherr << "Swath " << datasetname.c_str () << " cannot be attached.";
136  throw InternalErr (__FILE__, __LINE__, eherr.str ());
137  }
138 
139 
140  int32 tmp_rank = -1;
141  vector<int32>tmp_dims;
142  tmp_dims.resize(rank);
143  char tmp_dimlist[1024];
144  int32 type =-1;
145  intn r = -1;
146  r = fieldinfofunc (swathid, const_cast < char *>(fieldname.c_str ()),
147  &tmp_rank, tmp_dims.data(), &type, tmp_dimlist);
148  if (r != 0) {
149  detachfunc (swathid);
150  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
151  ostringstream eherr;
152  eherr << "Field " << fieldname.c_str () << " information cannot be obtained.";
153  throw InternalErr (__FILE__, __LINE__, eherr.str ());
154  }
155 
156  vector<int32> newdims;
157  newdims.resize(tmp_rank);
158 
159 
160  switch (type) {
161  case DFNT_INT8:
162  {
163  if(true == no_interpolation) {
164  vector<int8>val;
165  val.resize(nelms);
166  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
167  if (r != 0) {
168  detachfunc (swathid);
169  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
170  ostringstream eherr;
171  eherr << "field " << fieldname.c_str () << "cannot be read.";
172  throw InternalErr (__FILE__, __LINE__, eherr.str ());
173  }
174 #ifndef SIGNED_BYTE_TO_INT32
175  set_value ((dods_byte *) val.data(), nelms);
176 #else
177  vector<int32>newval;
178  newval.resize(nelms);
179 
180  for (int counter = 0; counter < nelms; counter++)
181  newval[counter] = (int32) (val[counter]);
182  set_value ((dods_int32 *) newval.data(), nelms);
183 #endif
184  }
185  else {
186  // Obtaining the total value and interpolating the data
187  // according to dimension map
188  vector <int8> total_val8;
189  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val8, newdims);
190 
191  if (r != 0) {
192  detachfunc (swathid);
193  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
194  ostringstream eherr;
195  eherr << "field " << fieldname.c_str () << "cannot be read.";
196  throw InternalErr (__FILE__, __LINE__, eherr.str ());
197 
198  }
199  vector<int8>val8;
200  val8.resize(nelms);
201  Field2DSubset (val8.data(), newdims[0], newdims[1],total_val8.data(),
202  offset32.data(), count32.data(), step32.data());
203 
204 #ifndef SIGNED_BYTE_TO_INT32
205  set_value ((dods_byte *) val8.data(), nelms);
206 #else
207  vector<int32>newval;
208  newval.resize(nelms);
209 
210  for (int counter = 0; counter < nelms; counter++)
211  newval[counter] = (int32) (val8[counter]);
212  set_value ((dods_int32 *) newval.data(), nelms);
213 #endif
214  }
215  }
216  break;
217  case DFNT_UINT8:
218  case DFNT_UCHAR8:
219  {
220  if(no_interpolation == false) {
221  vector <uint8> total_val_uint8;
222  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_uint8, newdims);
223 
224  if (r != 0) {
225  detachfunc (swathid);
226  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
227  ostringstream eherr;
228  eherr << "field " << fieldname.c_str () << "cannot be read.";
229  throw InternalErr (__FILE__, __LINE__, eherr.str ());
230 
231  }
232  vector<uint8>val_uint8;
233  val_uint8.resize(nelms);
234  Field2DSubset (val_uint8.data(), newdims[0], newdims[1],total_val_uint8.data(),
235  offset32.data(), count32.data(), step32.data());
236 
237  set_value ((dods_byte *) val_uint8.data(), nelms);
238  }
239  else {
240 
241  vector<uint8>val;
242  val.resize(nelms);
243  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
244  if (r != 0) {
245  detachfunc (swathid);
246  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
247  ostringstream eherr;
248  eherr << "field " << fieldname.c_str () << "cannot be read.";
249  throw InternalErr (__FILE__, __LINE__, eherr.str ());
250  }
251 
252  set_value ((dods_byte *) val.data(), nelms);
253  }
254  }
255  break;
256 
257  case DFNT_INT16:
258  {
259  if(no_interpolation == false) {
260  vector <int16> total_val_int16;
261  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_int16, newdims);
262 
263  if (r != 0) {
264  detachfunc (swathid);
265  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
266  ostringstream eherr;
267  eherr << "field " << fieldname.c_str () << "cannot be read.";
268  throw InternalErr (__FILE__, __LINE__, eherr.str ());
269 
270  }
271  vector<int16>val_int16;
272  val_int16.resize(nelms);
273  Field2DSubset (val_int16.data(), newdims[0], newdims[1],total_val_int16.data(),
274  offset32.data(), count32.data(), step32.data());
275 
276  set_value ((dods_int16 *) val_int16.data(), nelms);
277  }
278 
279  else {
280  vector<int16>val;
281  val.resize(nelms);
282  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
283  if (r != 0) {
284  detachfunc (swathid);
285  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
286  ostringstream eherr;
287  eherr << "field " << fieldname.c_str () << "cannot be read.";
288  throw InternalErr (__FILE__, __LINE__, eherr.str ());
289  }
290  set_value ((dods_int16 *) val.data(), nelms);
291  }
292  }
293  break;
294 
295  case DFNT_UINT16:
296  {
297  if(no_interpolation == false) {
298  vector <uint16> total_val_uint16;
299  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_uint16, newdims);
300 
301  if (r != 0) {
302  detachfunc (swathid);
303  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
304  ostringstream eherr;
305  eherr << "field " << fieldname.c_str () << "cannot be read.";
306  throw InternalErr (__FILE__, __LINE__, eherr.str ());
307 
308  }
309  vector<uint16>val_uint16;
310  val_uint16.resize(nelms);
311  Field2DSubset (val_uint16.data(), newdims[0], newdims[1],total_val_uint16.data(),
312  offset32.data(), count32.data(), step32.data());
313 
314  set_value ((dods_uint16 *) val_uint16.data(), nelms);
315  }
316 
317  else {
318  vector<uint16>val;
319  val.resize(nelms);
320  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
321  if (r != 0) {
322  detachfunc (swathid);
323  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
324  ostringstream eherr;
325  eherr << "field " << fieldname.c_str () << "cannot be read.";
326  throw InternalErr (__FILE__, __LINE__, eherr.str ());
327  }
328 
329  set_value ((dods_uint16 *) val.data(), nelms);
330  }
331  }
332  break;
333 
334  case DFNT_INT32:
335  {
336  if(no_interpolation == false) {
337  vector <int32> total_val_int32;
338  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_int32, newdims);
339 
340  if (r != 0) {
341  detachfunc (swathid);
342  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
343  ostringstream eherr;
344  eherr << "field " << fieldname.c_str () << "cannot be read.";
345  throw InternalErr (__FILE__, __LINE__, eherr.str ());
346 
347  }
348  vector<int32>val_int32;
349  val_int32.resize(nelms);
350  Field2DSubset (val_int32.data(), newdims[0], newdims[1],total_val_int32.data(),
351  offset32.data(), count32.data(), step32.data());
352 
353  set_value ((dods_int32 *) val_int32.data(), nelms);
354  }
355  else {
356  vector<int32>val;
357  val.resize(nelms);
358  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
359  if (r != 0) {
360  detachfunc (swathid);
361  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
362  ostringstream eherr;
363  eherr << "field " << fieldname.c_str () << "cannot be read.";
364  throw InternalErr (__FILE__, __LINE__, eherr.str ());
365  }
366 
367  set_value ((dods_int32 *) val.data(), nelms);
368  }
369  }
370  break;
371 
372  case DFNT_UINT32:
373  {
374  if(no_interpolation == false) {
375  vector <uint32> total_val_uint32;
376  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_uint32, newdims);
377 
378  if (r != 0) {
379  detachfunc (swathid);
380  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
381  ostringstream eherr;
382  eherr << "field " << fieldname.c_str () << "cannot be read.";
383  throw InternalErr (__FILE__, __LINE__, eherr.str ());
384 
385  }
386  vector<uint32>val_uint32;
387  val_uint32.resize(nelms);
388  Field2DSubset (val_uint32.data(), newdims[0], newdims[1],total_val_uint32.data(),
389  offset32.data(), count32.data(), step32.data());
390 
391  set_value ((dods_uint32 *) val_uint32.data(), nelms);
392  }
393  else {
394  vector<uint32>val;
395  val.resize(nelms);
396  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
397  if (r != 0) {
398  detachfunc (swathid);
399  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
400  ostringstream eherr;
401  eherr << "field " << fieldname.c_str () << "cannot be read.";
402  throw InternalErr (__FILE__, __LINE__, eherr.str ());
403  }
404 
405  set_value ((dods_uint32 *) val.data(), nelms);
406  }
407  }
408  break;
409  case DFNT_FLOAT32:
410  {
411  if(no_interpolation == false) {
412  vector <float32> total_val_f32;
413  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_f32, newdims);
414 
415  if (r != 0) {
416  detachfunc (swathid);
417  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
418  ostringstream eherr;
419  eherr << "field " << fieldname.c_str () << "cannot be read.";
420  throw InternalErr (__FILE__, __LINE__, eherr.str ());
421 
422  }
423  vector<float32>val_f32;
424  val_f32.resize(nelms);
425  Field2DSubset (val_f32.data(), newdims[0], newdims[1],total_val_f32.data(),
426  offset32.data(), count32.data(), step32.data());
427 
428  set_value ((dods_float32 *) val_f32.data(), nelms);
429  }
430  else {
431  vector<float32>val;
432  val.resize(nelms);
433  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
434  if (r != 0) {
435  detachfunc (swathid);
436  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
437  ostringstream eherr;
438  eherr << "field " << fieldname.c_str () << "cannot be read.";
439  throw InternalErr (__FILE__, __LINE__, eherr.str ());
440  }
441  set_value ((dods_float32 *) val.data(), nelms);
442  }
443  }
444  break;
445  case DFNT_FLOAT64:
446  {
447  if(no_interpolation == false) {
448  vector <float64> total_val_f64;
449  r = GetFieldValue (swathid, fieldname, dm_dimsizes,dm_offsets,dm_incs, total_val_f64, newdims);
450 
451  if (r != 0) {
452  detachfunc (swathid);
453  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
454  ostringstream eherr;
455  eherr << "field " << fieldname.c_str () << "cannot be read.";
456  throw InternalErr (__FILE__, __LINE__, eherr.str ());
457 
458  }
459  vector<float64>val_f64;
460  val_f64.resize(nelms);
461  Field2DSubset (val_f64.data(), newdims[0], newdims[1],total_val_f64.data(),
462  offset32.data(), count32.data(), step32.data());
463 
464  set_value ((dods_float64 *) val_f64.data(), nelms);
465  }
466  else {
467  vector<float64>val;
468  val.resize(nelms);
469  r = readfieldfunc (swathid, const_cast < char *>(fieldname.c_str ()), offset32.data(), step32.data(), count32.data(), val.data());
470  if (r != 0) {
471  detachfunc (swathid);
472  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
473  ostringstream eherr;
474  eherr << "field " << fieldname.c_str () << "cannot be read.";
475  throw InternalErr (__FILE__, __LINE__, eherr.str ());
476  }
477 
478  set_value ((dods_float64 *) val.data(), nelms);
479  }
480  }
481  break;
482  default:
483  detachfunc (swathid);
484  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
485  throw InternalErr (__FILE__, __LINE__, "unsupported data type.");
486  }
487 
488  r = detachfunc (swathid);
489  if (r != 0) {
490  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
491  ostringstream eherr;
492  eherr << "Swath " << datasetname.c_str () << " cannot be detached.";
493  throw InternalErr (__FILE__, __LINE__, eherr.str ());
494  }
495 
496  HDFCFUtil::close_fileid(-1,-1,-1,sfid,check_pass_fileid_key);
497 
498 #if 0
499  r = closefunc (sfid);
500  if (r != 0) {
501  ostringstream eherr;
502  eherr << "Swath " << filename.c_str () << " cannot be closed.";
503  throw InternalErr (__FILE__, __LINE__, eherr.str ());
504  }
505 #endif
506 
507  return false;
508 }
509 
510 // Standard way of DAP handlers to pass the coordinates of the subsetted region to the handlers
511 // Return the number of elements to read.
512 int
513 HDFEOS2ArraySwathGeoMultiDimMapField::format_constraint (int *offset, int *step, int *count)
514 {
515  int nels = 1;
516  int id = 0;
517 
518  Dim_iter p = dim_begin ();
519  while (p != dim_end ()) {
520 
521  int start = dimension_start (p, true);
522  int stride = dimension_stride (p, true);
523  int stop = dimension_stop (p, true);
524 
525  // Check for illegal constraint
526  if (start > stop) {
527  ostringstream oss;
528  oss << "Array/Grid hyperslab start point "<< start <<
529  " is greater than stop point " << stop <<".";
530  throw Error(malformed_expr, oss.str());
531  }
532 
533  offset[id] = start;
534  step[id] = stride;
535  count[id] = ((stop - start) / stride) + 1; // count of elements
536  nels *= count[id]; // total number of values for variable
537 
538  BESDEBUG ("h4",
539  "=format_constraint():"
540  << "id=" << id << " offset=" << offset[id]
541  << " step=" << step[id]
542  << " count=" << count[id]
543  << endl);
544 
545  id++;
546  p++;
547  }// end of while
548 
549  return nels;
550 }
551 
552 // Subset of latitude and longitude to follow the parameters
553 // from the DAP expression constraint
554 template < class T >
555 bool HDFEOS2ArraySwathGeoMultiDimMapField::Field2DSubset (T * outlatlon,
556  const int ,
557  const int minordim,
558  T * latlon,
559  const int32 * offset,
560  const int32 * count,
561  const int32 * step) const
562 {
563 #if 0
564  T (*templatlonptr)[majordim][minordim] = (T *[][]) latlon;
565 #endif
566  int i = 0;
567  int j = 0;
568 
569  // do subsetting
570  // Find the correct index
571  int dim0count = count[0];
572  int dim1count = count[1];
573 
574  int dim0index[dim0count];
575  int dim1index[dim1count];
576 
577  for (i = 0; i < count[0]; i++) // count[0] is the least changing dimension
578  dim0index[i] = offset[0] + i * step[0];
579 
580 
581  for (j = 0; j < count[1]; j++)
582  dim1index[j] = offset[1] + j * step[1];
583 
584  // Now assign the subsetting data
585  int k = 0;
586 
587  for (i = 0; i < count[0]; i++) {
588  for (j = 0; j < count[1]; j++) {
589 #if 0
590  outlatlon[k] = (*templatlonptr)[dim0index[i]][dim1index[j]];
591 #endif
592  outlatlon[k] = *(latlon + (dim0index[i] * minordim) + dim1index[j]);
593  k++;
594  }
595  }
596  return true;
597 }
598 
599 // Get latitude and longitude fields.
600 // It will call expand_dimmap_field to interpolate latitude and longitude.
601 template < class T > int
602 HDFEOS2ArraySwathGeoMultiDimMapField::
603 GetFieldValue (int32 swathid, const string & geofieldname,
604  const vector <int>&dimsizes,const vector<int>&offset,const vector<int>&inc,
605  vector < T > &vals, vector<int32>&newdims)
606 {
607 
608  int32 ret = -1;
609  int32 size = -1;
610  int32 sw_rank = -1;
611  int32 dims[130];
612  int32 type = -1;
613 
614  // Two dimensions for lat/lon; each dimension name is < 64 characters,
615  // The dimension names are separated by a comma.
616  char dimlist[130];
617  ret = SWfieldinfo (swathid, const_cast < char *>(geofieldname.c_str ()),
618  &sw_rank, dims, &type, dimlist);
619  if (ret != 0)
620  return -1;
621 
622  if(sw_rank !=2)
623  return -1;
624 
625  size = 1;
626  for (int i = 0; i <sw_rank; i++)
627  size *= dims[i];
628 
629  vals.resize (size);
630 
631  ret = SWreadfield (swathid, const_cast < char *>(geofieldname.c_str ()),
632  NULL, NULL, NULL, (void *) vals.data());
633  if (ret != 0)
634  return -1;
635 
636  vector < string > dimname;
637  HDFCFUtil::Split (dimlist, ',', dimname);
638 
639  for (int i = 0; i < sw_rank; i++) {
640 
641 #if 0
642 cerr<<"dimsizes["<<i<<"]: " <<dimsizes[i]<<endl;
643 cerr<<"offset["<<i<<"]: " <<offset[i]<<endl;
644 cerr<<"inc["<<i<<"]: " <<inc[i]<<endl;
645 #endif
646  int r;
647  r = _expand_dimmap_field (&vals, sw_rank, dims, i, dimsizes[i], offset[i], inc[i]);
648  if (r != 0)
649  return -1;
650  }
651 
652  // dims[] are expanded already.
653  for (int i = 0; i < sw_rank; i++) {
654  //cerr<<"i "<< i << " "<< dims[i] <<endl;
655  if (dims[i] < 0)
656  return -1;
657  newdims[i] = dims[i];
658  }
659 
660  return 0;
661 }
662 
663 // expand the dimension map field.
664 template < class T > int
665 HDFEOS2ArraySwathGeoMultiDimMapField::_expand_dimmap_field (vector < T >
666  *pvals, int32 sw_rank,
667  int32 dimsa[],
668  int dimindex,
669  int32 ddimsize,
670  int32 offset,
671  int32 inc) const
672 {
673  vector < T > orig = *pvals;
674  vector < int32 > pos;
675  vector < int32 > dims;
676  vector < int32 > newdims;
677  pos.resize (sw_rank);
678  dims.resize (sw_rank);
679 
680  for (int i = 0; i < sw_rank; i++) {
681  pos[i] = 0;
682  dims[i] = dimsa[i];
683  }
684  newdims = dims;
685  newdims[dimindex] = ddimsize;
686  dimsa[dimindex] = ddimsize;
687 
688  int newsize = 1;
689 
690  for (int i = 0; i < sw_rank; i++) {
691  newsize *= newdims[i];
692  }
693  pvals->clear ();
694  pvals->resize (newsize);
695 
696  for (;;) {
697  // if end
698  if (pos[0] == dims[0]) {
699  // we past then end
700  break;
701  }
702  else if (pos[dimindex] == 0) {
703  // extract 1D values
704  vector < T > v;
705  for (int i = 0; i < dims[dimindex]; i++) {
706  pos[dimindex] = i;
707  v.push_back (orig[INDEX_nD_TO_1D (dims, pos)]);
708  }
709  // expand them
710 
711  vector < T > w;
712  for (int32 j = 0; j < ddimsize; j++) {
713  int32 i = (j - offset) / inc;
714  T f;
715 
716  if (i * inc + offset == j) // perfect match
717  {
718  f = (v[i]);
719  }
720  else {
721  int32 i1 = 0;
722  int32 i2 = (i<=0)?1:0;
723  int32 j1 = 0;
724  int32 j2 = 0;
725 
726 #if 0
727  if (i <= 0) {
728  //i1 = 0;
729  i2 = 1;
730  }
731 #endif
732  if ((unsigned int) i + 1 >= v.size ()) {
733  i1 = v.size () - 2;
734  i2 = v.size () - 1;
735  }
736  else {
737  i1 = i;
738  i2 = i + 1;
739  }
740  j1 = i1 * inc + offset;
741  j2 = i2 * inc + offset;
742  f = (((j - j1) * v[i2] + (j2 - j) * v[i1]) / (j2 - j1));
743  }
744  w.push_back (f);
745  pos[dimindex] = j;
746  (*pvals)[INDEX_nD_TO_1D (newdims, pos)] = f;
747  }
748  pos[dimindex] = 0;
749  }
750  // next pos
751  pos[sw_rank - 1]++;
752  for (int i = sw_rank - 1; i > 0; i--) {
753  if (pos[i] == dims[i]) {
754  pos[i] = 0;
755  pos[i - 1]++;
756  }
757  }
758  }
759 
760  return 0;
761 }
762 
763 
764 #endif
HDFCFUtil::Split
static void Split(const char *s, int len, char sep, std::vector< std::string > &names)
Definition: HDFCFUtil.cc:82
HDFCFUtil::close_fileid
static void close_fileid(int32 sdfd, int32 file_id, int32 gridfd, int32 swathfd, bool pass_fileid_key)
Definition: HDFCFUtil.cc:3650