bes  Updated for version 3.20.10
HDFCFUtil.cc
1 #include "config.h"
2 #include "config_hdf.h"
3 
4 #include "HDFCFUtil.h"
5 #include <BESDebug.h>
6 #include <BESLog.h>
7 #include <math.h>
8 #include"dodsutil.h"
9 #include"HDFSPArray_RealField.h"
10 #include"HDFSPArrayMissField.h"
11 #include"HDFEOS2GeoCFProj.h"
12 #include"HDFEOS2GeoCF1D.h"
13 
14 #include <libdap/debug.h>
15 
16 #define SIGNED_BYTE_TO_INT32 1
17 
18 // HDF datatype headers for both the default and the CF options
19 #include "HDFByte.h"
20 #include "HDFInt16.h"
21 #include "HDFUInt16.h"
22 #include "HDFInt32.h"
23 #include "HDFUInt32.h"
24 #include "HDFFloat32.h"
25 #include "HDFFloat64.h"
26 #include "HDFStr.h"
27 #include "HDF4RequestHandler.h"
28 //
29 using namespace std;
30 using namespace libdap;
31 
32 #define ERR_LOC1(x) #x
33 #define ERR_LOC2(x) ERR_LOC1(x)
34 #define ERR_LOC __FILE__ " : " ERR_LOC2(__LINE__)
35 
36 // Check the BES key.
37 // This function will check a BES key specified at the file h4.conf.in.
38 // If the key's value is either true or yes. The handler claims to find
39 // a key and will do some operations. Otherwise, will do different operations.
40 // For example, One may find a line H4.EnableCF=true at h4.conf.in.
41 // That means, the HDF4 handler will handle the HDF4 files by following CF conventions.
42 #if 0
43 bool
44 HDFCFUtil::check_beskeys(const string& key) {
45 
46  bool found = false;
47  string doset ="";
48  const string dosettrue ="true";
49  const string dosetyes = "yes";
50 
51  TheBESKeys::TheKeys()->get_value( key, doset, found ) ;
52  if( true == found ) {
53  doset = BESUtil::lowercase( doset ) ;
54  if( dosettrue == doset || dosetyes == doset )
55  return true;
56  }
57  return false;
58 
59 }
60 #endif
61 
66 static void
67 split_helper(vector<string> &tokens, const string &text, const char sep)
68 {
69  string::size_type start = 0, end = 0;
70  while ((end = text.find(sep, start)) != string::npos) {
71  tokens.push_back(text.substr(start, end - start));
72  start = end + 1;
73  }
74  tokens.push_back(text.substr(start));
75 }
76 
77 // From a string separated by a separator to a list of string,
78 // for example, split "ab,c" to {"ab","c"}
79 void
80 HDFCFUtil::Split(const char *s, int len, char sep, std::vector<std::string> &names)
81 {
82  names.clear();
83  split_helper(names, string(s, len), sep);
84 #if 0
85  // Replaced with the above since valgrind reports errors
86  // with this code. jhrg 6/22/15
87  for (int i = 0, j = 0; j <= len; ++j) {
88  if ((j == len && len) || s[j] == sep) {
89  string elem(s + i, j - i);
90  names.push_back(elem);
91  i = j + 1;
92  continue;
93  }
94  }
95 #endif
96 }
97 
98 // Assume sz is Null terminated string.
99 void
100 HDFCFUtil::Split(const char *sz, char sep, std::vector<std::string> &names)
101 {
102  names.clear();
103 
104  // Split() was showing up in some valgrind runs as having an off-by-one
105  // error. I added this help track it down.
106  DBG(std::cerr << "HDFCFUtil::Split: sz: <" << sz << ">, sep: <" << sep << ">" << std::endl);
107 
108  split_helper(names, string(sz), sep);
109 #if 0
110  // Replaced with a direct call to the new helper code.
111  // jhrg 6/22/15
112  Split(sz, (int)strlen(sz), sep, names);
113 #endif
114 }
115 
116 #if 0
117 // From a string separated by a separator to a list of string,
118 // for example, split "ab,c" to {"ab","c"}
119 void
120 HDFCFUtil::Split(const char *s, int len, char sep, std::vector<std::string> &names)
121 {
122  names.clear();
123  for (int i = 0, j = 0; j <= len; ++j) {
124  if ((j == len && len) || s[j] == sep) {
125  string elem(s + i, j - i);
126  names.push_back(elem);
127  i = j + 1;
128  continue;
129  }
130  }
131 }
132 
133 // Assume sz is Null terminated string.
134 void
135 HDFCFUtil::Split(const char *sz, char sep, std::vector<std::string> &names)
136 {
137  Split(sz, (int)strlen(sz), sep, names);
138 }
139 
140 #endif
141 
142 // This is a safer way to insert and update a c++ map value.
143 // Otherwise, the local testsuite at The HDF Group will fail for HDF-EOS2 data
144 // under iMac machine platform.
145 // The implementation replaces the element even if the key exists.
146 // This function is equivalent to map[key]=value
147 bool
148 HDFCFUtil::insert_map(std::map<std::string,std::string>& m, string key, string val)
149 {
150  pair<map<string,string>::iterator, bool> ret;
151  ret = m.insert(make_pair(key, val));
152  if(ret.second == false){
153  m.erase(key);
154  ret = m.insert(make_pair(key, val));
155  if(ret.second == false){
156  BESDEBUG("h4","insert_map():insertion failed on Key=" << key << " Val=" << val << endl);
157  }
158  }
159  return ret.second;
160 }
161 
162 // Obtain CF string
163 string
165 {
166 
167  if(""==s)
168  return s;
169  string insertString(1,'_');
170 
171  // Always start with _ if the first character is not a letter
172  if (true == isdigit(s[0]))
173  s.insert(0,insertString);
174 
175  // New name conventions drop the first '/' from the path.
176  if ('/' ==s[0])
177  s.erase(0,1);
178 
179  for(unsigned int i=0; i < s.length(); i++)
180  if((false == isalnum(s[i])) && (s[i]!='_'))
181  s[i]='_';
182 
183  return s;
184 
185 }
186 
187 // Obtain the unique name for the clashed names and save it to set namelist.
188 // This is a recursive call. A unique name list is represented as a set.
189 // If we find that a name already exists in the nameset, we will add a number
190 // at the end of the name to form a new name. If the new name still exists
191 // in the nameset, we will increase the index number and check again until
192 // a unique name is generated.
193 void
194 HDFCFUtil::gen_unique_name(string &str,set<string>& namelist, int&clash_index) {
195 
196  pair<set<string>::iterator,bool> ret;
197  string newstr = "";
198  stringstream sclash_index;
199  sclash_index << clash_index;
200  newstr = str + sclash_index.str();
201 
202  ret = namelist.insert(newstr);
203  if (false == ret.second) {
204  clash_index++;
205  gen_unique_name(str,namelist,clash_index);
206  }
207  else
208  str = newstr;
209 }
210 
211 // General routine to handle the name clashing
212 // The input parameters include:
213 // name vector -> newobjnamelist(The name vector is changed to a unique name list
214 // a pre-allocated object name set ->objnameset(can be used to determine if a name exists)
215 void
216 HDFCFUtil::Handle_NameClashing(vector<string>&newobjnamelist,set<string>&objnameset) {
217 
218  pair<set<string>::iterator,bool> setret;
219  set<string>::iterator iss;
220 
221  vector<string> clashnamelist;
222  vector<string>::iterator ivs;
223 
224  // clash index to original index mapping
225  map<int,int> cl_to_ol;
226  int ol_index = 0;
227  int cl_index = 0;
228 
229  vector<string>::const_iterator irv;
230 
231  for (irv = newobjnamelist.begin(); irv != newobjnamelist.end(); ++irv) {
232  setret = objnameset.insert(*irv);
233  if (false == setret.second ) {
234  clashnamelist.insert(clashnamelist.end(),(*irv));
235  cl_to_ol[cl_index] = ol_index;
236  cl_index++;
237  }
238  ol_index++;
239  }
240 
241  // Now change the clashed elements to unique elements,
242  // Generate the set which has the same size as the original vector.
243  for (ivs=clashnamelist.begin(); ivs!=clashnamelist.end(); ivs++) {
244  int clash_index = 1;
245  string temp_clashname = *ivs +'_';
246  HDFCFUtil::gen_unique_name(temp_clashname,objnameset,clash_index);
247  *ivs = temp_clashname;
248  }
249 
250  // Now go back to the original vector, make it unique.
251  for (unsigned int i =0; i <clashnamelist.size(); i++)
252  newobjnamelist[cl_to_ol[i]] = clashnamelist[i];
253 
254 }
255 
256 // General routine to handle the name clashing
257 // The input parameter just includes:
258 // name vector -> newobjnamelist(The name vector is changed to a unique name list
259 void
260 HDFCFUtil::Handle_NameClashing(vector<string>&newobjnamelist) {
261 
262  set<string> objnameset;
263  Handle_NameClashing(newobjnamelist,objnameset);
264 }
265 
266 // Borrowed codes from ncdas.cc in netcdf_handle4 module.
267 string
268 HDFCFUtil::print_attr(int32 type, int loc, void *vals)
269 {
270  ostringstream rep;
271 
272  union {
273  char *cp;
274  unsigned char *ucp;
275  short *sp;
276  unsigned short *usp;
277  int32 /*nclong*/ *lp;
278  unsigned int *ui;
279  float *fp;
280  double *dp;
281  } gp;
282 
283  switch (type) {
284 
285  // Mapping both DFNT_UINT8 and DFNT_INT8 to unsigned char
286  // may cause overflow. Documented at jira ticket HFRHANDLER-169.
287  case DFNT_UINT8:
288  {
289  unsigned char uc;
290  gp.ucp = (unsigned char *) vals;
291 
292  uc = *(gp.ucp+loc);
293  rep << (int)uc;
294  return rep.str();
295  }
296 
297  case DFNT_INT8:
298  {
299  char c;
300  gp.cp = (char *) vals;
301 
302  c = *(gp.cp+loc);
303  rep << (int)c;
304  return rep.str();
305  }
306 
307  case DFNT_UCHAR:
308  case DFNT_CHAR:
309  {
310  // Use the customized escattr function. Don't escape \n,\t and \r. KY 2013-10-14
311  return escattr(static_cast<const char*>(vals));
312  }
313 
314  case DFNT_INT16:
315  {
316  gp.sp = (short *) vals;
317  rep << *(gp.sp+loc);
318  return rep.str();
319  }
320 
321  case DFNT_UINT16:
322  {
323  gp.usp = (unsigned short *) vals;
324  rep << *(gp.usp+loc);
325  return rep.str();
326  }
327 
328  case DFNT_INT32:
329  {
330  gp.lp = (int32 *) vals;
331  rep << *(gp.lp+loc);
332  return rep.str();
333  }
334 
335  case DFNT_UINT32:
336  {
337  gp.ui = (unsigned int *) vals;
338  rep << *(gp.ui+loc);
339  return rep.str();
340  }
341 
342  case DFNT_FLOAT:
343  {
344  float attr_val = *(float*)vals;
345  bool is_a_fin = isfinite(attr_val);
346  gp.fp = (float *) vals;
347  rep << showpoint;
348  // setprecision seeme to cause the one-bit error when
349  // converting from float to string. Watch whether this
350  // is an isue.
351  rep << setprecision(10);
352  rep << *(gp.fp+loc);
353  string tmp_rep_str = rep.str();
354  if (tmp_rep_str.find('.') == string::npos
355  && tmp_rep_str.find('e') == string::npos
356  && tmp_rep_str.find('E') == string::npos) {
357  if(true == is_a_fin)
358  rep << ".";
359  }
360  return rep.str();
361  }
362 
363  case DFNT_DOUBLE:
364  {
365 
366  double attr_val = *(double*)vals;
367  bool is_a_fin = isfinite(attr_val);
368  gp.dp = (double *) vals;
369  rep << std::showpoint;
370  rep << std::setprecision(17);
371  rep << *(gp.dp+loc);
372  string tmp_rep_str = rep.str();
373  if (tmp_rep_str.find('.') == string::npos
374  && tmp_rep_str.find('e') == string::npos
375  && tmp_rep_str.find('E') == string::npos) {
376  if(true == is_a_fin)
377  rep << ".";
378  }
379  return rep.str();
380  }
381  default:
382  return string("UNKNOWN");
383  }
384 
385 }
386 
387 // Print datatype in string. This is used to generate DAS.
388 string
390 {
391 
392  // I expanded the list based on libdap/AttrTable.h.
393  static const char UNKNOWN[]="Unknown";
394  static const char BYTE[]="Byte";
395  static const char INT16[]="Int16";
396  static const char UINT16[]="UInt16";
397  static const char INT32[]="Int32";
398  static const char UINT32[]="UInt32";
399  static const char FLOAT32[]="Float32";
400  static const char FLOAT64[]="Float64";
401  static const char STRING[]="String";
402 
403  // I got different cases from hntdefs.h.
404  switch (type) {
405 
406  case DFNT_CHAR:
407  return STRING;
408 
409  case DFNT_UCHAR8:
410  return STRING;
411 
412  case DFNT_UINT8:
413  return BYTE;
414 
415  case DFNT_INT8:
416 // ADD the macro
417  {
418 #ifndef SIGNED_BYTE_TO_INT32
419  return BYTE;
420 #else
421  return INT32;
422 #endif
423  }
424  case DFNT_UINT16:
425  return UINT16;
426 
427  case DFNT_INT16:
428  return INT16;
429 
430  case DFNT_INT32:
431  return INT32;
432 
433  case DFNT_UINT32:
434  return UINT32;
435 
436  case DFNT_FLOAT:
437  return FLOAT32;
438 
439  case DFNT_DOUBLE:
440  return FLOAT64;
441 
442  default:
443  return UNKNOWN;
444  }
445 
446 }
447 
448 // Obtain HDF4 datatype size.
449 short
451 {
452 
453  // .
454  switch (type) {
455 
456  case DFNT_CHAR:
457  return sizeof(char);
458 
459  case DFNT_UCHAR8:
460  return sizeof(unsigned char);
461 
462  case DFNT_UINT8:
463  return sizeof(unsigned char);
464 
465  case DFNT_INT8:
466 // ADD the macro
467  {
468 #ifndef SIGNED_BYTE_TO_INT32
469  return sizeof(char);
470 #else
471  return sizeof(int);
472 #endif
473  }
474  case DFNT_UINT16:
475  return sizeof(unsigned short);
476 
477  case DFNT_INT16:
478  return sizeof(short);
479 
480  case DFNT_INT32:
481  return sizeof(int);
482 
483  case DFNT_UINT32:
484  return sizeof(unsigned int);
485 
486  case DFNT_FLOAT:
487  return sizeof(float);
488 
489  case DFNT_DOUBLE:
490  return sizeof(double);
491 
492  default:
493  return -1;
494  }
495 
496 }
497 // Subset of latitude and longitude to follow the parameters from the DAP expression constraint
498 // Somehow this function doesn't work. Now it is not used. Still keep it here for the future investigation.
499 template < typename T >
500 void HDFCFUtil::LatLon2DSubset (T * outlatlon,
501  int majordim,
502  int minordim,
503  T * latlon,
504  int32 * offset,
505  int32 * count,
506  int32 * step)
507 {
508 
509  // float64 templatlon[majordim][minordim];
510 #if 0
511  // --std=c++11 on OSX causes 'typeof' to fail. This is a GNU gcc-specific
512  // keyword. jhrg 3/28/19
513  T (*templatlonptr)[majordim][minordim] = (typeof templatlonptr) latlon;
514 #endif
515  T (*templatlonptr)[majordim][minordim] = (T *) latlon;
516  int i, j, k;
517 
518  // do subsetting
519  // Find the correct index
520  int dim0count = count[0];
521  int dim1count = count[1];
522  int dim0index[dim0count], dim1index[dim1count];
523 
524  for (i = 0; i < count[0]; i++) // count[0] is the least changing dimension
525  dim0index[i] = offset[0] + i * step[0];
526 
527 
528  for (j = 0; j < count[1]; j++)
529  dim1index[j] = offset[1] + j * step[1];
530 
531  // Now assign the subsetting data
532  k = 0;
533 
534  for (i = 0; i < count[0]; i++) {
535  for (j = 0; j < count[1]; j++) {
536  outlatlon[k] = (*templatlonptr)[dim0index[i]][dim1index[j]];
537  k++;
538 
539  }
540  }
541 }
542 
543 // CF requires the _FillValue attribute datatype is the same as the corresponding field datatype.
544 // For some NASA files, this is not true.
545 // So we need to check if the _FillValue's datatype is the same as the attribute's.
546 // If not, we need to correct them.
547 void HDFCFUtil::correct_fvalue_type(AttrTable *at,int32 dtype) {
548 
549  AttrTable::Attr_iter it = at->attr_begin();
550  bool find_fvalue = false;
551  while (it!=at->attr_end() && false==find_fvalue) {
552  if (at->get_name(it) =="_FillValue")
553  {
554  find_fvalue = true;
555  string fillvalue ="";
556  string fillvalue_type ="";
557  if((*at->get_attr_vector(it)).size() !=1)
558  throw InternalErr(__FILE__,__LINE__,"The number of _FillValue must be 1.");
559  fillvalue = (*at->get_attr_vector(it)->begin());
560  fillvalue_type = at->get_type(it);
561  string var_type = HDFCFUtil::print_type(dtype);
562 
563  if(fillvalue_type != var_type){
564 
565  at->del_attr("_FillValue");
566 
567  if (fillvalue_type == "String") {
568 
569  // String fillvalue is always represented as /+octal numbers when its type is forced to
570  // change to string(check HFRHANDLER-223). So we have to retrieve it back.
571  if(fillvalue.size() >1) {
572 
573  long int fillvalue_int = 0;
574  vector<char> fillvalue_temp(fillvalue.size());
575  char *pEnd;
576  fillvalue_int = strtol((fillvalue.substr(1)).c_str(),&pEnd,8);
577  stringstream convert_str;
578  convert_str << fillvalue_int;
579  at->append_attr("_FillValue",var_type,convert_str.str());
580  }
581  else {
582 
583  // If somehow the fillvalue type is DFNT_CHAR or DFNT_UCHAR, and the size is 1,
584  // that means the fillvalue type is wrongly defined, we treat as a 8-bit integer number.
585  // Note, the code will only assume the value ranges from 0 to 128.(JIRA HFRHANDLER-248).
586  // KY 2014-04-24
587 
588  short fillvalue_int = fillvalue.at(0);
589 
590  stringstream convert_str;
591  convert_str << fillvalue_int;
592  if(fillvalue_int <0 || fillvalue_int >128)
593  throw InternalErr(__FILE__,__LINE__,
594  "If the fillvalue is a char type, the value must be between 0 and 128.");
595 
596 
597  at->append_attr("_FillValue",var_type,convert_str.str());
598  }
599  }
600 
601  else
602  at->append_attr("_FillValue",var_type,fillvalue);
603  }
604  }
605  it++;
606  }
607 
608 }
609 
610 // CF requires the scale_factor and add_offset attribute datatypes must be the same as the corresponding field datatype.
611 // For some NASA files, this is not true.
612 // So we need to check if the _FillValue's datatype is the same as the attribute's.
613 // If not, we need to correct them.
615 
616  AttrTable::Attr_iter it = at->attr_begin();
617  bool find_scale = false;
618  bool find_offset = false;
619 
620  // Declare scale_factor,add_offset, fillvalue and valid_range type in string format.
621  string scale_factor_type;
622  string add_offset_type;
623  string scale_factor_value;
624  string add_offset_value;
625 
626  while (it!=at->attr_end() &&((find_scale!=true) ||(find_offset!=true))) {
627  if (at->get_name(it) =="scale_factor")
628  {
629  find_scale = true;
630  scale_factor_value = (*at->get_attr_vector(it)->begin());
631  scale_factor_type = at->get_type(it);
632  }
633 
634  if(at->get_name(it)=="add_offset")
635  {
636  find_offset = true;
637  add_offset_value = (*at->get_attr_vector(it)->begin());
638  add_offset_type = at->get_type(it);
639  }
640 
641  it++;
642  }
643 
644  // Change offset type to the scale type
645  if((true==find_scale) && (true==find_offset)) {
646  if(scale_factor_type != add_offset_type) {
647  at->del_attr("add_offset");
648  at->append_attr("add_offset",scale_factor_type,add_offset_value);
649  }
650  }
651 }
652 
653 
654 #ifdef USE_HDFEOS2_LIB
655 
656 // For MODIS (confirmed by level 1B) products, values between 65500(MIN_NON_SCALE_SPECIAL_VALUE)
657 // and 65535(MAX_NON_SCALE_SPECIAL_VALUE) are treated as
658 // special values. These values represent non-physical data values caused by various failures.
659 // For example, 65533 represents "when Detector is saturated".
660 bool HDFCFUtil::is_special_value(int32 dtype, float fillvalue, float realvalue) {
661 
662  bool ret_value = false;
663 
664  if (DFNT_UINT16 == dtype) {
665 
666  int fillvalue_int = (int)fillvalue;
667  if (MAX_NON_SCALE_SPECIAL_VALUE == fillvalue_int) {
668  int realvalue_int = (int)realvalue;
669  if (realvalue_int <= MAX_NON_SCALE_SPECIAL_VALUE && realvalue_int >=MIN_NON_SCALE_SPECIAL_VALUE)
670  ret_value = true;
671  }
672  }
673 
674  return ret_value;
675 
676 }
677 
678 // Check if the MODIS file has dimension map and return the number of dimension maps
679 // Note: This routine is only applied to a MODIS geo-location file when the corresponding
680 // MODIS swath uses dimension maps and has the MODIS geo-location file under the same
681 // file directory. This is also restricted by turning on H4.EnableCheckMODISGeoFile to be true(file h4.conf.in).
682 // By default, this key is turned off. Also this function is only used once in one service. So it won't
683 // affect performance. KY 2014-02-18
684 int HDFCFUtil::check_geofile_dimmap(const string & geofilename) {
685 
686  int32 fileid = SWopen(const_cast<char*>(geofilename.c_str()),DFACC_READ);
687  if (fileid < 0)
688  return -1;
689  string swathname = "MODIS_Swath_Type_GEO";
690  int32 datasetid = SWattach(fileid,const_cast<char*>(swathname.c_str()));
691  if (datasetid < 0) {
692  SWclose(fileid);
693  return -1;
694  }
695 
696  int32 nummaps = 0;
697  int32 bufsize;
698 
699  // Obtain number of dimension maps and the buffer size.
700  if ((nummaps = SWnentries(datasetid, HDFE_NENTMAP, &bufsize)) == -1) {
701  SWdetach(datasetid);
702  SWclose(fileid);
703  return -1;
704  }
705 
706  SWdetach(datasetid);
707  SWclose(fileid);
708  return nummaps;
709 
710 }
711 
712 // Check if we need to change the datatype for MODIS fields. The datatype needs to be changed
713 // mainly because of non-CF scale and offset rules. To avoid violating CF conventions, we apply
714 // the non-CF MODIS scale and offset rule to MODIS data. So the final data type may be different
715 // than the original one due to this operation. For example, the original datatype may be int16.
716 // After applying the scale/offset rule, the datatype may become float32.
717 // The following are useful notes about MODIS SCALE OFFSET HANDLING
718 // Note: MODIS Scale and offset handling needs to re-organized. But it may take big efforts.
719 // Instead, I remove the global variable mtype, and _das; move the old calculate_dtype code
720 // back to HDFEOS2.cc. The code is a little better organized. If possible, we may think to overhaul
721 // the handling of MODIS scale-offset part. KY 2012-6-19
722 //
723 bool HDFCFUtil::change_data_type(DAS & das, SOType scaletype, const string &new_field_name)
724 {
725 
726  AttrTable *at = das.get_table(new_field_name);
727 
728  // The following codes do these:
729  // For MODIS level 1B(using the swath name), check radiance,reflectance etc.
730  // If the DISABLESCALE key is true, no need to check scale and offset for other types.
731  // Otherwise, continue checking the scale and offset names.
732  // KY 2013-12-13
733 
734  if(scaletype!=DEFAULT_CF_EQU && at!=NULL)
735  {
736  AttrTable::Attr_iter it = at->attr_begin();
737  string scale_factor_value="";
738  string add_offset_value="0";
739  string radiance_scales_value="";
740  string radiance_offsets_value="";
741  string reflectance_scales_value="";
742  string reflectance_offsets_value="";
743  string scale_factor_type, add_offset_type;
744  while (it!=at->attr_end())
745  {
746  if(at->get_name(it)=="radiance_scales")
747  radiance_scales_value = *(at->get_attr_vector(it)->begin());
748  if(at->get_name(it)=="radiance_offsets")
749  radiance_offsets_value = *(at->get_attr_vector(it)->begin());
750  if(at->get_name(it)=="reflectance_scales")
751  reflectance_scales_value = *(at->get_attr_vector(it)->begin());
752  if(at->get_name(it)=="reflectance_offsets")
753  reflectance_offsets_value = *(at->get_attr_vector(it)->begin());
754 
755  // Ideally may just check if the attribute name is scale_factor.
756  // But don't know if some products use attribut name like "scale_factor "
757  // So now just filter out the attribute scale_factor_err. KY 2012-09-20
758  if(at->get_name(it).find("scale_factor")!=string::npos){
759  string temp_attr_name = at->get_name(it);
760  if (temp_attr_name != "scale_factor_err") {
761  scale_factor_value = *(at->get_attr_vector(it)->begin());
762  scale_factor_type = at->get_type(it);
763  }
764  }
765  if(at->get_name(it).find("add_offset")!=string::npos)
766  {
767  string temp_attr_name = at->get_name(it);
768  if (temp_attr_name !="add_offset_err") {
769  add_offset_value = *(at->get_attr_vector(it)->begin());
770  add_offset_type = at->get_type(it);
771  }
772  }
773  it++;
774  }
775 
776  if((radiance_scales_value.length()!=0 && radiance_offsets_value.length()!=0)
777  || (reflectance_scales_value.length()!=0 && reflectance_offsets_value.length()!=0))
778  return true;
779 
780  if(scale_factor_value.length()!=0)
781  {
782  if(!(atof(scale_factor_value.c_str())==1 && atof(add_offset_value.c_str())==0))
783  return true;
784  }
785  }
786 
787  return false;
788 }
789 
790 // Obtain the MODIS swath dimension map info.
791 void HDFCFUtil::obtain_dimmap_info(const string& filename,HDFEOS2::Dataset*dataset,
792  vector<struct dimmap_entry> & dimmaps,
793  string & modis_geofilename, bool& geofile_has_dimmap) {
794 
795 
796  HDFEOS2::SwathDataset *sw = static_cast<HDFEOS2::SwathDataset *>(dataset);
797  const vector<HDFEOS2::SwathDataset::DimensionMap*>& origdimmaps = sw->getDimensionMaps();
798  vector<HDFEOS2::SwathDataset::DimensionMap*>::const_iterator it_dmap;
799  struct dimmap_entry tempdimmap;
800 
801  // if having dimension maps, we need to retrieve the dimension map info.
802  for(size_t i=0;i<origdimmaps.size();i++){
803  tempdimmap.geodim = origdimmaps[i]->getGeoDimension();
804  tempdimmap.datadim = origdimmaps[i]->getDataDimension();
805  tempdimmap.offset = origdimmaps[i]->getOffset();
806  tempdimmap.inc = origdimmaps[i]->getIncrement();
807  dimmaps.push_back(tempdimmap);
808  }
809 
810 #if 0
811  string check_modis_geofile_key ="H4.EnableCheckMODISGeoFile";
812  bool check_geofile_key = false;
813  check_geofile_key = HDFCFUtil::check_beskeys(check_modis_geofile_key);
814 #endif
815 
816  // Only when there is dimension map, we need to consider the additional MODIS geolocation files.
817  // Will check if the check modis_geo_location file key is turned on.
818  if((origdimmaps.size() != 0) && (true == HDF4RequestHandler::get_enable_check_modis_geo_file()) ) {
819 
820  // Has to use C-style since basename and dirname are not C++ routines.
821  char*tempcstr;
822  tempcstr = new char [filename.size()+1];
823  strncpy (tempcstr,filename.c_str(),filename.size());
824  string basefilename = basename(tempcstr);
825  string dirfilename = dirname(tempcstr);
826  delete [] tempcstr;
827 
828  // If the current file is a MOD03 or a MYD03 file, we don't need to check the extra MODIS geolocation file at all.
829  bool is_modis_geofile = false;
830  if(basefilename.size() >5) {
831  if((0 == basefilename.compare(0,5,"MOD03")) || (0 == basefilename.compare(0,5,"MYD03")))
832  is_modis_geofile = true;
833  }
834 
835  // This part is implemented specifically for supporting MODIS dimension map data.
836  // MODIS Aqua Swath dimension map geolocation file always starts with MYD03
837  // MODIS Terra Swath dimension map geolocation file always starts with MOD03
838  // We will check the first three characters to see if the dimension map geolocation file exists.
839  // An example MODIS swath filename is MOD05_L2.A2008120.0000.005.2008121182723.hdf
840  // An example MODIS geo-location file name is MOD03.A2008120.0000.005.2010003235220.hdf
841  // Notice that the "A2008120.0000" in the middle of the name is the "Acquistion Date" of the data
842  // So the geo-location file name should have exactly the same string. We will use this
843  // string to identify if a MODIS geo-location file exists or not.
844  // Note the string size is 14(.A2008120.0000).
845  // More information of naming convention, check http://modis-atmos.gsfc.nasa.gov/products_filename.html
846  // KY 2010-5-10
847 
848 
849  // Obtain string "MYD" or "MOD"
850  // Here we need to consider when MOD03 or MYD03 use the dimension map.
851  if ((false == is_modis_geofile) && (basefilename.size() >3)) {
852 
853  string fnameprefix = basefilename.substr(0,3);
854 
855  if(fnameprefix == "MYD" || fnameprefix =="MOD") {
856  size_t fnamemidpos = basefilename.find(".A");
857  if(fnamemidpos != string::npos) {
858  string fnamemiddle = basefilename.substr(fnamemidpos,14);
859  if(fnamemiddle.size()==14) {
860  string geofnameprefix = fnameprefix+"03";
861  // geofnamefp will be something like "MOD03.A2008120.0000"
862  string geofnamefp = geofnameprefix + fnamemiddle;
863  DIR *dirp;
864  struct dirent* dirs;
865 
866  // Go through the directory to see if we have the corresponding MODIS geolocation file
867  dirp = opendir(dirfilename.c_str());
868  if (NULL == dirp)
869  throw InternalErr(__FILE__,__LINE__,"opendir fails.");
870 
871  while ((dirs = readdir(dirp))!= NULL){
872  if(strncmp(dirs->d_name,geofnamefp.c_str(),geofnamefp.size())==0){
873  modis_geofilename = dirfilename + "/"+ dirs->d_name;
874  int num_dimmap = HDFCFUtil::check_geofile_dimmap(modis_geofilename);
875  if (num_dimmap < 0) {
876  closedir(dirp);
877  throw InternalErr(__FILE__,__LINE__,"this file is not a MODIS geolocation file.");
878  }
879  geofile_has_dimmap = (num_dimmap >0)?true:false;
880  break;
881  }
882  }
883  closedir(dirp);
884  }
885  }
886  }
887  }
888  }
889 }
890 
891 // This is for the case that the separate MODIS geo-location file is used.
892 // Some geolocation names at the MODIS data file are not consistent with
893 // the names in the MODIS geo-location file. So need to correct them.
894 bool HDFCFUtil::is_modis_dimmap_nonll_field(string & fieldname) {
895 
896  bool modis_dimmap_nonll_field = false;
897  vector<string> modis_dimmap_nonll_fieldlist;
898 
899  modis_dimmap_nonll_fieldlist.push_back("Height");
900  modis_dimmap_nonll_fieldlist.push_back("SensorZenith");
901  modis_dimmap_nonll_fieldlist.push_back("SensorAzimuth");
902  modis_dimmap_nonll_fieldlist.push_back("Range");
903  modis_dimmap_nonll_fieldlist.push_back("SolarZenith");
904  modis_dimmap_nonll_fieldlist.push_back("SolarAzimuth");
905  modis_dimmap_nonll_fieldlist.push_back("Land/SeaMask");
906  modis_dimmap_nonll_fieldlist.push_back("gflags");
907  modis_dimmap_nonll_fieldlist.push_back("Solar_Zenith");
908  modis_dimmap_nonll_fieldlist.push_back("Solar_Azimuth");
909  modis_dimmap_nonll_fieldlist.push_back("Sensor_Azimuth");
910  modis_dimmap_nonll_fieldlist.push_back("Sensor_Zenith");
911 
912  map<string,string>modis_field_to_geofile_field;
913  map<string,string>::iterator itmap;
914  modis_field_to_geofile_field["Solar_Zenith"] = "SolarZenith";
915  modis_field_to_geofile_field["Solar_Azimuth"] = "SolarAzimuth";
916  modis_field_to_geofile_field["Sensor_Zenith"] = "SensorZenith";
917  modis_field_to_geofile_field["Solar_Azimuth"] = "SolarAzimuth";
918 
919  for (unsigned int i = 0; i <modis_dimmap_nonll_fieldlist.size(); i++) {
920 
921  if (fieldname == modis_dimmap_nonll_fieldlist[i]) {
922  itmap = modis_field_to_geofile_field.find(fieldname);
923  if (itmap !=modis_field_to_geofile_field.end())
924  fieldname = itmap->second;
925  modis_dimmap_nonll_field = true;
926  break;
927  }
928  }
929 
930  return modis_dimmap_nonll_field;
931 }
932 
933 void HDFCFUtil::add_scale_offset_attrs(AttrTable*at,const std::string& s_type, float svalue_f, double svalue_d, bool add_offset_found,
934  const std::string& o_type, float ovalue_f, double ovalue_d) {
935  at->del_attr("scale_factor");
936  string print_rep;
937  if(s_type!="Float64") {
938  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&svalue_f));
939  at->append_attr("scale_factor", "Float32", print_rep);
940  }
941  else {
942  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT64,0,(void*)(&svalue_d));
943  at->append_attr("scale_factor", "Float64", print_rep);
944  }
945 
946  if (true == add_offset_found) {
947  at->del_attr("add_offset");
948  if(o_type!="Float64") {
949  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&ovalue_f));
950  at->append_attr("add_offset", "Float32", print_rep);
951  }
952  else {
953  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT64,0,(void*)(&ovalue_d));
954  at->append_attr("add_offset", "Float64", print_rep);
955  }
956  }
957 }
958 
959 void HDFCFUtil::add_scale_str_offset_attrs(AttrTable*at,const std::string& s_type, const std::string& s_value_str, bool add_offset_found,
960  const std::string& o_type, float ovalue_f, double ovalue_d) {
961  at->del_attr("scale_factor");
962  string print_rep;
963  if(s_type!="Float64")
964  at->append_attr("scale_factor", "Float32", s_value_str);
965  else
966  at->append_attr("scale_factor", "Float64", s_value_str);
967 
968  if (true == add_offset_found) {
969  at->del_attr("add_offset");
970  if(o_type!="Float64") {
971  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&ovalue_f));
972  at->append_attr("add_offset", "Float32", print_rep);
973  }
974  else {
975  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT64,0,(void*)(&ovalue_d));
976  at->append_attr("add_offset", "Float64", print_rep);
977  }
978  }
979 }
981 void HDFCFUtil::handle_modis_special_attrs_disable_scale_comp(AttrTable *at,
982  const string &filename,
983  bool is_grid,
984  const string & newfname,
985  SOType sotype){
986 
987  // Declare scale_factor,add_offset, fillvalue and valid_range type in string format.
988  string scale_factor_type;
989  string add_offset_type;
990 
991  // Scale and offset values
992  string scale_factor_value="";
993  float orig_scale_value_float = 1;
994  float orig_scale_value_double = 1;
995  string add_offset_value="0";
996  float orig_offset_value_float = 0;
997  float orig_offset_value_double = 0;
998  bool add_offset_found = false;
999 
1000 
1001  // Go through all attributes to find scale_factor,add_offset,_FillValue,valid_range
1002  // Number_Type information
1003  AttrTable::Attr_iter it = at->attr_begin();
1004  while (it!=at->attr_end())
1005  {
1006  if(at->get_name(it)=="scale_factor")
1007  {
1008  scale_factor_value = (*at->get_attr_vector(it)->begin());
1009  scale_factor_type = at->get_type(it);
1010  if(scale_factor_type =="Float64")
1011  orig_scale_value_double=atof(scale_factor_value.c_str());
1012  else
1013  orig_scale_value_float = atof(scale_factor_value.c_str());
1014  }
1015 
1016  if(at->get_name(it)=="add_offset")
1017  {
1018  add_offset_value = (*at->get_attr_vector(it)->begin());
1019  add_offset_type = at->get_type(it);
1020 
1021  if(add_offset_type == "Float64")
1022  orig_offset_value_double = atof(add_offset_value.c_str());
1023  else
1024  orig_offset_value_float = atof(add_offset_value.c_str());
1025  add_offset_found = true;
1026  }
1027 
1028  it++;
1029  }
1030 
1031  // According to our observations, it seems that MODIS products ALWAYS use the "scale" factor to
1032  // make bigger values smaller.
1033  // So for MODIS_MUL_SCALE products, if the scale of some variable is greater than 1,
1034  // it means that for this variable, the MODIS type for this variable may be MODIS_DIV_SCALE.
1035  // For the similar logic, we may need to change MODIS_DIV_SCALE to MODIS_MUL_SCALE and MODIS_EQ_SCALE
1036  // to MODIS_DIV_SCALE.
1037  // We indeed find such a case. HDF-EOS2 Grid MODIS_Grid_1km_2D of MOD(or MYD)09GA is a MODIS_EQ_SCALE.
1038  // However,
1039  // the scale_factor of the variable Range_1 in the MOD09GA product is 25. According to our observation,
1040  // this variable should be MODIS_DIV_SCALE.We verify this is true according to MODIS 09 product document
1041  // http://modis-sr.ltdri.org/products/MOD09_UserGuide_v1_3.pdf.
1042  // Since this conclusion is based on our observation, we would like to add a BESlog to detect if we find
1043  // the similar cases so that we can verify with the corresponding product documents to see if this is true.
1044  // More information,
1045  // We just verified with the data producer, the scale_factor for Range_1 and Range_c is 25 but the
1046  // equation is still multiplication instead of division. So we have to make this as a special case that
1047  // we don't want to change the scale and offset settings.
1048  // KY 2014-01-13
1049 
1050 
1051  if(scale_factor_value.length()!=0) {
1052  if (MODIS_EQ_SCALE == sotype || MODIS_MUL_SCALE == sotype) {
1053  if (orig_scale_value_float > 1 || orig_scale_value_double >1) {
1054 
1055  bool need_change_scale = true;
1056  if(true == is_grid) {
1057  if ((filename.size() >5) && ((filename.compare(0,5,"MOD09") == 0)|| (filename.compare(0,5,"MYD09")==0))) {
1058  if ((newfname.size() >5) && newfname.find("Range") != string::npos)
1059  need_change_scale = false;
1060  }
1061  else if((filename.size() >7)&&((filename.compare(0,7,"MOD16A2") == 0)|| (filename.compare(0,7,"MYD16A2")==0) ||
1062  (filename.compare(0,7,"MOD16A3")==0) || (filename.compare(0,7,"MYD16A3")==0)))
1063  need_change_scale = false;
1064  }
1065  if(true == need_change_scale) {
1066  sotype = MODIS_DIV_SCALE;
1067  (*BESLog::TheLog())<< "The field " << newfname << " scale factor is "<< scale_factor_value << endl
1068  << " But the original scale factor type is MODIS_MUL_SCALE or MODIS_EQ_SCALE. " << endl
1069  << " Now change it to MODIS_DIV_SCALE. "<<endl;
1070  }
1071  }
1072  }
1073 
1074  if (MODIS_DIV_SCALE == sotype) {
1075  if (orig_scale_value_float < 1 || orig_scale_value_double<1) {
1076  sotype = MODIS_MUL_SCALE;
1077  (*BESLog::TheLog())<< "The field " << newfname << " scale factor is "<< scale_factor_value << endl
1078  << " But the original scale factor type is MODIS_DIV_SCALE. " << endl
1079  << " Now change it to MODIS_MUL_SCALE. "<<endl;
1080  }
1081  }
1082 
1083 
1084  if ((MODIS_MUL_SCALE == sotype) &&(true == add_offset_found)) {
1085 
1086  float new_offset_value_float=0;
1087  double new_offset_value_double=0;
1088  if(add_offset_type!="Float64")
1089  new_offset_value_float = (orig_offset_value_float ==0)?0:(-1 * orig_offset_value_float *orig_scale_value_float);
1090  else
1091  new_offset_value_double = (orig_offset_value_double ==0)?0:(-1 * orig_offset_value_double *orig_scale_value_double);
1092 
1093  // May need to use another function to avoid the rounding error by atof.
1094  add_scale_str_offset_attrs(at,scale_factor_type,scale_factor_value,add_offset_found,
1095  add_offset_type,new_offset_value_float,new_offset_value_double);
1096 
1097  }
1098 
1099  if (MODIS_DIV_SCALE == sotype) {
1100 
1101  float new_scale_value_float=1;
1102  double new_scale_value_double=1;
1103  float new_offset_value_float=0;
1104  double new_offset_value_double=0;
1105 
1106 
1107  if(scale_factor_type !="Float64") {
1108  new_scale_value_float = 1.0/orig_scale_value_float;
1109  if (true == add_offset_found) {
1110  if(add_offset_type !="Float64")
1111  new_offset_value_float = (orig_offset_value_float==0)?0:(-1 * orig_offset_value_float *new_scale_value_float);
1112  else
1113  new_offset_value_double = (orig_offset_value_double==0)?0:(-1 * orig_offset_value_double *new_scale_value_float);
1114  }
1115  }
1116 
1117  else {
1118  new_scale_value_double = 1.0/orig_scale_value_double;
1119  if (true == add_offset_found) {
1120  if(add_offset_type !="Float64")
1121  new_offset_value_float = (orig_offset_value_float==0)?0:(-1 * orig_offset_value_float *new_scale_value_double);
1122  else
1123  new_offset_value_double = (orig_offset_value_double==0)?0:(-1 * orig_offset_value_double *new_scale_value_double);
1124  }
1125  }
1126 
1127  add_scale_offset_attrs(at,scale_factor_type,new_scale_value_float,new_scale_value_double,add_offset_found,
1128  add_offset_type,new_offset_value_float,new_offset_value_double);
1129 
1130  }
1131 
1132  }
1133 
1134 }
1135 
1136 // These routines will handle scale_factor,add_offset,valid_min,valid_max and other attributes
1137 // such as Number_Type to make sure the CF is followed.
1138 // For example, For the case that the scale and offset rule doesn't follow CF, the scale_factor and add_offset attributes are renamed
1139 // to orig_scale_factor and orig_add_offset to keep the original field info.
1140 // Note: This routine should only be called when MODIS Scale and offset rules don't follow CF.
1141 // Input parameters:
1142 // AttrTable at - DAS attribute table
1143 // string newfname - field name: this is used to identify special MODIS level 1B emissive and refSB fields
1144 // SOType sotype - MODIS scale and offset type
1145 // bool gridname_change_valid_range - Flag to check if this is the special MODIS VIP product
1146 // bool changedtype - Flag to check if the datatype of this field needs to be changed
1147 // bool change_fvtype - Flag to check if the attribute _FillValue needs to change to data type
1148 
1150 // Divide this function into smaller functions:
1151 //
1152 void HDFCFUtil::handle_modis_special_attrs(AttrTable *at, const string & filename,
1153  bool is_grid,const string & newfname,
1154  SOType sotype, bool gridname_change_valid_range,
1155  bool changedtype, bool & change_fvtype){
1156 
1157  // Declare scale_factor,add_offset, fillvalue and valid_range type in string format.
1158  string scale_factor_type;
1159  string add_offset_type;
1160  string fillvalue_type;
1161  string valid_range_type;
1162 
1163 
1164  // Scale and offset values
1165  string scale_factor_value="";
1166  float orig_scale_value = 1;
1167  string add_offset_value="0";
1168  float orig_offset_value = 0;
1169  bool add_offset_found = false;
1170 
1171  // Fillvalue
1172  string fillvalue="";
1173 
1174  // Valid range value
1175  string valid_range_value="";
1176 
1177  bool has_valid_range = false;
1178 
1179  // We may need to change valid_range to valid_min and valid_max. Initialize them here.
1180  float orig_valid_min = 0;
1181  float orig_valid_max = 0;
1182 
1183  // Number_Type also needs to be adjusted when datatype is changed
1184  string number_type_value="";
1185  string number_type_dap_type="";
1186 
1187  // Go through all attributes to find scale_factor,add_offset,_FillValue,valid_range
1188  // Number_Type information
1189  AttrTable::Attr_iter it = at->attr_begin();
1190  while (it!=at->attr_end())
1191  {
1192  if(at->get_name(it)=="scale_factor")
1193  {
1194  scale_factor_value = (*at->get_attr_vector(it)->begin());
1195  orig_scale_value = atof(scale_factor_value.c_str());
1196  scale_factor_type = at->get_type(it);
1197  }
1198 
1199  if(at->get_name(it)=="add_offset")
1200  {
1201  add_offset_value = (*at->get_attr_vector(it)->begin());
1202  orig_offset_value = atof(add_offset_value.c_str());
1203  add_offset_type = at->get_type(it);
1204  add_offset_found = true;
1205  }
1206 
1207  if(at->get_name(it)=="_FillValue")
1208  {
1209  fillvalue = (*at->get_attr_vector(it)->begin());
1210  fillvalue_type = at->get_type(it);
1211  }
1212 
1213  if(at->get_name(it)=="valid_range")
1214  {
1215  vector<string> *avalue = at->get_attr_vector(it);
1216  vector<string>::iterator ait = avalue->begin();
1217  while(ait!=avalue->end())
1218  {
1219  valid_range_value += *ait;
1220  ait++;
1221  if(ait!=avalue->end())
1222  valid_range_value += ", ";
1223  }
1224  valid_range_type = at->get_type(it);
1225  if (false == gridname_change_valid_range) {
1226  orig_valid_min = (float)(atof((avalue->at(0)).c_str()));
1227  orig_valid_max = (float)(atof((avalue->at(1)).c_str()));
1228  }
1229  has_valid_range = true;
1230  }
1231 
1232  if(true == changedtype && (at->get_name(it)=="Number_Type"))
1233  {
1234  number_type_value = (*at->get_attr_vector(it)->begin());
1235  number_type_dap_type= at->get_type(it);
1236  }
1237 
1238  it++;
1239  }
1240 
1241  // Rename scale_factor and add_offset attribute names. Otherwise, they will be
1242  // misused by CF tools to generate wrong data values based on the CF scale and offset rule.
1243  if(scale_factor_value.length()!=0)
1244  {
1245  if(!(atof(scale_factor_value.c_str())==1 && atof(add_offset_value.c_str())==0)) //Rename them.
1246  {
1247  at->del_attr("scale_factor");
1248  at->append_attr("orig_scale_factor", scale_factor_type, scale_factor_value);
1249  if(add_offset_found)
1250  {
1251  at->del_attr("add_offset");
1252  at->append_attr("orig_add_offset", add_offset_type, add_offset_value);
1253  }
1254  }
1255  }
1256 
1257  // Change _FillValue datatype
1258  if(true == changedtype && fillvalue.length()!=0 && fillvalue_type!="Float32" && fillvalue_type!="Float64")
1259  {
1260  change_fvtype = true;
1261  at->del_attr("_FillValue");
1262  at->append_attr("_FillValue", "Float32", fillvalue);
1263  }
1264 
1265  float valid_max = 0;
1266  float valid_min = 0;
1267 
1268  it = at->attr_begin();
1269  bool handle_modis_l1b = false;
1270 
1271  // MODIS level 1B's Emissive and RefSB fields' scale_factor and add_offset
1272  // get changed according to different vertical levels.
1273  // So we need to handle them specifically.
1274  if (sotype == MODIS_MUL_SCALE && true ==changedtype) {
1275 
1276  // Emissive should be at the end of the field name such as "..._Emissive"
1277  string emissive_str = "Emissive";
1278  string RefSB_str = "RefSB";
1279  bool is_emissive_field = false;
1280  bool is_refsb_field = false;
1281  if(newfname.find(emissive_str)!=string::npos) {
1282  if(0 == newfname.compare(newfname.size()-emissive_str.size(),emissive_str.size(),emissive_str))
1283  is_emissive_field = true;
1284  }
1285 
1286  if(newfname.find(RefSB_str)!=string::npos) {
1287  if(0 == newfname.compare(newfname.size()-RefSB_str.size(),RefSB_str.size(),RefSB_str))
1288  is_refsb_field = true;
1289  }
1290 
1291  // We will calculate the maximum and minimum scales and offsets within all the vertical levels.
1292  if ((true == is_emissive_field) || (true== is_refsb_field)){
1293 
1294  float scale_max = 0;
1295  float scale_min = 100000;
1296 
1297  float offset_max = 0;
1298  float offset_min = 0;
1299 
1300  float temp_var_val = 0;
1301 
1302  string orig_long_name_value;
1303  string modify_long_name_value;
1304  string str_removed_from_long_name=" Scaled Integers";
1305  string radiance_units_value;
1306 
1307  while (it!=at->attr_end())
1308  {
1309  // Radiance field(Emissive field)
1310  if(true == is_emissive_field) {
1311 
1312  if ("radiance_scales" == (at->get_name(it))) {
1313  vector<string> *avalue = at->get_attr_vector(it);
1314  for (vector<string>::const_iterator ait = avalue->begin();ait !=avalue->end();++ait) {
1315  temp_var_val = (float)(atof((*ait).c_str()));
1316  if (temp_var_val > scale_max)
1317  scale_max = temp_var_val;
1318  if (temp_var_val < scale_min)
1319  scale_min = temp_var_val;
1320  }
1321  }
1322 
1323  if ("radiance_offsets" == (at->get_name(it))) {
1324  vector<string> *avalue = at->get_attr_vector(it);
1325  for (vector<string>::const_iterator ait = avalue->begin();ait !=avalue->end();++ait) {
1326  temp_var_val = (float)(atof((*ait).c_str()));
1327  if (temp_var_val > offset_max)
1328  offset_max = temp_var_val;
1329  if (temp_var_val < scale_min)
1330  offset_min = temp_var_val;
1331  }
1332  }
1333 
1334  if ("long_name" == (at->get_name(it))) {
1335  orig_long_name_value = (*at->get_attr_vector(it)->begin());
1336  if (orig_long_name_value.find(str_removed_from_long_name)!=string::npos) {
1337  if(0 == orig_long_name_value.compare(orig_long_name_value.size()-str_removed_from_long_name.size(),
1338  str_removed_from_long_name.size(),str_removed_from_long_name)) {
1339 
1340  modify_long_name_value =
1341  orig_long_name_value.substr(0,orig_long_name_value.size()-str_removed_from_long_name.size());
1342  at->del_attr("long_name");
1343  at->append_attr("long_name","String",modify_long_name_value);
1344  at->append_attr("orig_long_name","String",orig_long_name_value);
1345  }
1346  }
1347  }
1348 
1349  if ("radiance_units" == (at->get_name(it)))
1350  radiance_units_value = (*at->get_attr_vector(it)->begin());
1351  }
1352 
1353  if (true == is_refsb_field) {
1354  if ("reflectance_scales" == (at->get_name(it))) {
1355  vector<string> *avalue = at->get_attr_vector(it);
1356  for (vector<string>::const_iterator ait = avalue->begin();ait !=avalue->end();++ait) {
1357  temp_var_val = (float)(atof((*ait).c_str()));
1358  if (temp_var_val > scale_max)
1359  scale_max = temp_var_val;
1360  if (temp_var_val < scale_min)
1361  scale_min = temp_var_val;
1362  }
1363  }
1364 
1365  if ("reflectance_offsets" == (at->get_name(it))) {
1366 
1367  vector<string> *avalue = at->get_attr_vector(it);
1368  for (vector<string>::const_iterator ait = avalue->begin();ait !=avalue->end();++ait) {
1369  temp_var_val = (float)(atof((*ait).c_str()));
1370  if (temp_var_val > offset_max)
1371  offset_max = temp_var_val;
1372  if (temp_var_val < scale_min)
1373  offset_min = temp_var_val;
1374  }
1375  }
1376 
1377  if ("long_name" == (at->get_name(it))) {
1378  orig_long_name_value = (*at->get_attr_vector(it)->begin());
1379  if (orig_long_name_value.find(str_removed_from_long_name)!=string::npos) {
1380  if(0 == orig_long_name_value.compare(orig_long_name_value.size()-str_removed_from_long_name.size(),
1381  str_removed_from_long_name.size(),str_removed_from_long_name)) {
1382 
1383  modify_long_name_value =
1384  orig_long_name_value.substr(0,orig_long_name_value.size()-str_removed_from_long_name.size());
1385  at->del_attr("long_name");
1386  at->append_attr("long_name","String",modify_long_name_value);
1387  at->append_attr("orig_long_name","String",orig_long_name_value);
1388  }
1389  }
1390  }
1391  }
1392  it++;
1393  }
1394 
1395  // Calculate the final valid_max and valid_min.
1396  if (scale_min <= 0)
1397  throw InternalErr(__FILE__,__LINE__,"the scale factor should always be greater than 0.");
1398 
1399  if (orig_valid_max > offset_min)
1400  valid_max = (orig_valid_max-offset_min)*scale_max;
1401  else
1402  valid_max = (orig_valid_max-offset_min)*scale_min;
1403 
1404  if (orig_valid_min > offset_max)
1405  valid_min = (orig_valid_min-offset_max)*scale_min;
1406  else
1407  valid_min = (orig_valid_min -offset_max)*scale_max;
1408 
1409  // These physical variables should be greater than 0
1410  if (valid_min < 0)
1411  valid_min = 0;
1412 
1413  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&valid_min));
1414  at->append_attr("valid_min","Float32",print_rep);
1415  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&valid_max));
1416  at->append_attr("valid_max","Float32",print_rep);
1417  at->del_attr("valid_range");
1418  handle_modis_l1b = true;
1419 
1420  // Change the units for the emissive field
1421  if (true == is_emissive_field && radiance_units_value.size() >0) {
1422  at->del_attr("units");
1423  at->append_attr("units","String",radiance_units_value);
1424  }
1425  }
1426  }
1427 
1428  // Handle other valid_range attributes
1429  if(true == changedtype && true == has_valid_range && false == handle_modis_l1b) {
1430 
1431  // If the gridname_change_valid_range is true, call a special function to handle this.
1432  if (true == gridname_change_valid_range)
1433  HDFCFUtil::handle_modis_vip_special_attrs(valid_range_value,scale_factor_value,valid_min,valid_max);
1434  else if(scale_factor_value.length()!=0) {
1435 
1436  // We found MODIS products always scale to a smaller value. If somehow the original scale factor
1437  // is smaller than 1, the scale/offset should be the multiplication rule.(new_data =scale*(old_data-offset))
1438  // If the original scale factor is greater than 1, the scale/offset rule should be division rule.
1439  // New_data = (old_data-offset)/scale.
1440  // We indeed find such a case. HDF-EOS2 Grid MODIS_Grid_1km_2D of MOD(or MYD)09GA is a MODIS_EQ_SCALE.
1441  // However,
1442  // the scale_factor of the variable Range_1 in the MOD09GA product is 25. According to our observation,
1443  // this variable should be MODIS_DIV_SCALE.We verify this is true according to MODIS 09 product document
1444  // http://modis-sr.ltdri.org/products/MOD09_UserGuide_v1_3.pdf.
1445  // Since this conclusion is based on our observation, we would like to add a BESlog to detect if we find
1446  // the similar cases so that we can verify with the corresponding product documents to see if this is true.
1447  // More information,
1448  // We just verified with the data producer, the scale_factor for Range_1 and Range_c is 25 but the
1449  // equation is still multiplication instead of division. So we have to make this as a special case that
1450  // we don't want to change the scale and offset settings.
1451  // KY 2014-01-13
1452 
1453  if (MODIS_EQ_SCALE == sotype || MODIS_MUL_SCALE == sotype) {
1454  if (orig_scale_value > 1) {
1455 
1456  bool need_change_scale = true;
1457  if(true == is_grid) {
1458  if ((filename.size() >5) && ((filename.compare(0,5,"MOD09") == 0)|| (filename.compare(0,5,"MYD09")==0))) {
1459  if ((newfname.size() >5) && newfname.find("Range") != string::npos)
1460  need_change_scale = false;
1461  }
1462 
1463  else if((filename.size() >7)&&((filename.compare(0,7,"MOD16A2") == 0)|| (filename.compare(0,7,"MYD16A2")==0) ||
1464  (filename.compare(0,7,"MOD16A3")==0) || (filename.compare(0,7,"MYD16A3")==0)))
1465  need_change_scale = false;
1466  }
1467  if(true == need_change_scale) {
1468  sotype = MODIS_DIV_SCALE;
1469  (*BESLog::TheLog())<< "The field " << newfname << " scale factor is "<< orig_scale_value << endl
1470  << " But the original scale factor type is MODIS_MUL_SCALE or MODIS_EQ_SCALE. " << endl
1471  << " Now change it to MODIS_DIV_SCALE. "<<endl;
1472  }
1473  }
1474  }
1475 
1476  if (MODIS_DIV_SCALE == sotype) {
1477  if (orig_scale_value < 1) {
1478  sotype = MODIS_MUL_SCALE;
1479  (*BESLog::TheLog())<< "The field " << newfname << " scale factor is "<< orig_scale_value << endl
1480  << " But the original scale factor type is MODIS_DIV_SCALE. " << endl
1481  << " Now change it to MODIS_MUL_SCALE. "<<endl;
1482  }
1483  }
1484 
1485  if(sotype == MODIS_MUL_SCALE) {
1486  valid_min = (orig_valid_min - orig_offset_value)*orig_scale_value;
1487  valid_max = (orig_valid_max - orig_offset_value)*orig_scale_value;
1488  }
1489  else if (sotype == MODIS_DIV_SCALE) {
1490  valid_min = (orig_valid_min - orig_offset_value)/orig_scale_value;
1491  valid_max = (orig_valid_max - orig_offset_value)/orig_scale_value;
1492  }
1493  else if (sotype == MODIS_EQ_SCALE) {
1494  valid_min = orig_valid_min * orig_scale_value + orig_offset_value;
1495  valid_max = orig_valid_max * orig_scale_value + orig_offset_value;
1496  }
1497  }
1498  else {// This is our current assumption.
1499  if (MODIS_MUL_SCALE == sotype || MODIS_DIV_SCALE == sotype) {
1500  valid_min = orig_valid_min - orig_offset_value;
1501  valid_max = orig_valid_max - orig_offset_value;
1502  }
1503  else if (MODIS_EQ_SCALE == sotype) {
1504  valid_min = orig_valid_min + orig_offset_value;
1505  valid_max = orig_valid_max + orig_offset_value;
1506  }
1507  }
1508 
1509  // Append the corrected valid_min and valid_max. (valid_min,valid_max) is equivalent to valid_range.
1510  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&valid_min));
1511  at->append_attr("valid_min","Float32",print_rep);
1512  print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&valid_max));
1513  at->append_attr("valid_max","Float32",print_rep);
1514  at->del_attr("valid_range");
1515  }
1516 
1517  // We also find that there is an attribute called "Number_Type" that will stores the original attribute
1518  // datatype. If the datatype gets changed, the attribute is confusion. here we can change the attribute
1519  // name to "Number_Type_Orig"
1520  if(true == changedtype && number_type_dap_type !="" ) {
1521  at->del_attr("Number_Type");
1522  at->append_attr("Number_Type_Orig",number_type_dap_type,number_type_value);
1523  }
1524 }
1525 
1526  // This routine makes the MeaSUREs VIP attributes follow CF.
1527 // valid_range attribute uses char type instead of the raw data's datatype.
1528 void HDFCFUtil::handle_modis_vip_special_attrs(const std::string& valid_range_value,
1529  const std::string& scale_factor_value,
1530  float& valid_min, float &valid_max) {
1531 
1532  int16 vip_orig_valid_min = 0;
1533  int16 vip_orig_valid_max = 0;
1534 
1535  size_t found = valid_range_value.find_first_of(",");
1536  size_t found_from_end = valid_range_value.find_last_of(",");
1537  if (string::npos == found)
1538  throw InternalErr(__FILE__,__LINE__,"should find the separator ,");
1539  if (found != found_from_end)
1540  throw InternalErr(__FILE__,__LINE__,"There should be only one separator.");
1541 
1542 #if 0
1543  //istringstream(valid_range_value.substr(0,found))>>orig_valid_min;
1544  //istringstream(valid_range_value.substr(found+1))>>orig_valid_max;
1545 #endif
1546 
1547  vip_orig_valid_min = atoi((valid_range_value.substr(0,found)).c_str());
1548  vip_orig_valid_max = atoi((valid_range_value.substr(found+1)).c_str());
1549 
1550  int16 scale_factor_number = 1;
1551 
1552  scale_factor_number = atoi(scale_factor_value.c_str());
1553 
1554  if(scale_factor_number !=0) {
1555  valid_min = (float)(vip_orig_valid_min/scale_factor_number);
1556  valid_max = (float)(vip_orig_valid_max/scale_factor_number);
1557  }
1558  else
1559  throw InternalErr(__FILE__,__LINE__,"The scale_factor_number should not be zero.");
1560 }
1561 
1562 // Make AMSR-E attributes follow CF.
1563 // Change SCALE_FACTOR and OFFSET to CF names: scale_factor and add_offset.
1564 void HDFCFUtil::handle_amsr_attrs(AttrTable *at) {
1565 
1566  AttrTable::Attr_iter it = at->attr_begin();
1567  string scale_factor_value="", add_offset_value="0";
1568  string scale_factor_type, add_offset_type;
1569  bool OFFSET_found = false;
1570  bool Scale_found = false;
1571  bool SCALE_FACTOR_found = false;
1572 
1573  while (it!=at->attr_end())
1574  {
1575  if(at->get_name(it)=="SCALE_FACTOR")
1576  {
1577  scale_factor_value = (*at->get_attr_vector(it)->begin());
1578  scale_factor_type = at->get_type(it);
1579  SCALE_FACTOR_found = true;
1580  }
1581 
1582  if(at->get_name(it)=="Scale")
1583  {
1584  scale_factor_value = (*at->get_attr_vector(it)->begin());
1585  scale_factor_type = at->get_type(it);
1586  Scale_found = true;
1587  }
1588 
1589  if(at->get_name(it)=="OFFSET")
1590  {
1591  add_offset_value = (*at->get_attr_vector(it)->begin());
1592  add_offset_type = at->get_type(it);
1593  OFFSET_found = true;
1594  }
1595  it++;
1596  }
1597 
1598  if (true == SCALE_FACTOR_found) {
1599  at->del_attr("SCALE_FACTOR");
1600  at->append_attr("scale_factor",scale_factor_type,scale_factor_value);
1601  }
1602 
1603  if (true == Scale_found) {
1604  at->del_attr("Scale");
1605  at->append_attr("scale_factor",scale_factor_type,scale_factor_value);
1606  }
1607 
1608  if (true == OFFSET_found) {
1609  at->del_attr("OFFSET");
1610  at->append_attr("add_offset",add_offset_type,add_offset_value);
1611  }
1612 
1613 }
1614 
1615 //This function obtains the latitude and longitude dimension info. of an
1616 //HDF-EOS2 grid after the handler translates the HDF-EOS to CF.
1617 // Dimension info. includes dimension name and dimension size.
1618 void HDFCFUtil::obtain_grid_latlon_dim_info(HDFEOS2::GridDataset* gdset,
1619  string & dim0name,
1620  int32 & dim0size,
1621  string & dim1name,
1622  int32& dim1size){
1623 
1624  const vector<HDFEOS2::Field*>gfields = gdset->getDataFields();
1625  vector<HDFEOS2::Field*>::const_iterator it_gf;
1626  for (it_gf = gfields.begin();it_gf != gfields.end();++it_gf) {
1627 
1628  // Check the dimensions for Latitude
1629  if(1 == (*it_gf)->getFieldType()) {
1630  const vector<HDFEOS2::Dimension*>& dims= (*it_gf)->getCorrectedDimensions();
1631 
1632  //2-D latitude
1633  if(2 == dims.size()) {
1634  // Most time, it is YDim Major. We will see if we have a real case to
1635  // check if the handling is right for the XDim Major case.
1636  if(true == (*it_gf)->getYDimMajor()) {
1637  dim0name = dims[0]->getName();
1638  dim0size = dims[0]->getSize();
1639  dim1name = dims[1]->getName();
1640  dim1size = dims[1]->getSize();
1641  }
1642  else {
1643  dim0name = dims[1]->getName();
1644  dim0size = dims[1]->getSize();
1645  dim1name = dims[0]->getName();
1646  dim1size = dims[0]->getSize();
1647  }
1648  break;
1649  }
1650 
1651  //1-D latitude
1652  if( 1 == dims.size()) {
1653  dim0name = dims[0]->getName();
1654  dim0size = dims[0]->getSize();
1655  }
1656  }
1657 
1658  // Longitude, if longitude is checked first, it goes here.
1659  if(2 == (*it_gf)->getFieldType()) {
1660  const vector<HDFEOS2::Dimension*>& dims= (*it_gf)->getCorrectedDimensions();
1661  if(2 == dims.size()) {
1662 
1663  // Most time, it is YDim Major. We will see if we have a real case to
1664  // check if the handling is right for the XDim Major case.
1665  if(true == (*it_gf)->getYDimMajor()) {
1666  dim0name = dims[0]->getName();
1667  dim0size = dims[0]->getSize();
1668  dim1name = dims[1]->getName();
1669  dim1size = dims[1]->getSize();
1670  }
1671  else {
1672  dim0name = dims[1]->getName();
1673  dim0size = dims[1]->getSize();
1674  dim1name = dims[0]->getName();
1675  dim1size = dims[0]->getSize();
1676  }
1677  break;
1678  }
1679  if( 1 == dims.size()) {
1680  dim1name = dims[0]->getName();
1681  dim1size = dims[0]->getSize();
1682 
1683  }
1684  }
1685  }
1686  if(""==dim0name || ""==dim1name || dim0size<0 || dim1size <0)
1687  throw InternalErr (__FILE__, __LINE__,"Fail to obtain grid lat/lon dimension info.");
1688 
1689 }
1690 
1691 // This function adds the 1-D cf grid projection mapping attribute to data variables
1692 // it is called by the function add_cf_grid_attrs.
1693 void HDFCFUtil::add_cf_grid_mapping_attr(DAS &das, HDFEOS2::GridDataset*gdset,const string& cf_projection,
1694  const string & dim0name,int32 dim0size,const string &dim1name,int32 dim1size) {
1695 
1696  // Check >=2-D fields, check if they hold the dim0name,dim0size etc., yes, add the attribute cf_projection.
1697  const vector<HDFEOS2::Field*>gfields = gdset->getDataFields();
1698  vector<HDFEOS2::Field*>::const_iterator it_gf;
1699  for (it_gf = gfields.begin();it_gf != gfields.end();++it_gf) {
1700 
1701  if(0 == (*it_gf)->getFieldType() && (*it_gf)->getRank() >1) {
1702  bool has_dim0 = false;
1703  bool has_dim1 = false;
1704  const vector<HDFEOS2::Dimension*>& dims= (*it_gf)->getCorrectedDimensions();
1705  for (vector<HDFEOS2::Dimension *>::const_iterator j =
1706  dims.begin(); j!= dims.end();++j){
1707  if((*j)->getName()== dim0name && (*j)->getSize() == dim0size)
1708  has_dim0 = true;
1709  else if((*j)->getName()== dim1name && (*j)->getSize() == dim1size)
1710  has_dim1 = true;
1711 
1712  }
1713  if(true == has_dim0 && true == has_dim1) {// Need to add the grid_mapping attribute
1714  AttrTable *at = das.get_table((*it_gf)->getNewName());
1715  if (!at)
1716  at = das.add_table((*it_gf)->getNewName(), new AttrTable);
1717 
1718  // The dummy projection name is the value of the grid_mapping attribute
1719  at->append_attr("grid_mapping","String",cf_projection);
1720  }
1721  }
1722  }
1723 }
1724 
1725 //This function adds 1D grid mapping CF attributes to CV and data variables.
1726 void HDFCFUtil::add_cf_grid_cv_attrs(DAS & das, HDFEOS2::GridDataset *gdset) {
1727 
1728  //1. Check the projection information, now, we only handle sinusoidal now
1729  if(GCTP_SNSOID == gdset->getProjection().getCode()) {
1730 
1731  //2. Obtain the dimension information from latitude and longitude(fieldtype =1 or fieldtype =2)
1732  string dim0name,dim1name;
1733  int32 dim0size = -1,dim1size = -1;
1734  HDFCFUtil::obtain_grid_latlon_dim_info(gdset,dim0name,dim0size,dim1name,dim1size);
1735 
1736  //3. Add 1D CF attributes to the 1-D CV variables and the dummy projection variable
1737  AttrTable *at = das.get_table(dim0name);
1738  if (!at)
1739  at = das.add_table(dim0name, new AttrTable);
1740  at->append_attr("standard_name","String","projection_y_coordinate");
1741 
1742  string long_name="y coordinate of projection for grid "+ gdset->getName();
1743  at->append_attr("long_name","String",long_name);
1744  // Change this to meter.
1745  at->append_attr("units","string","meter");
1746 
1747  at->append_attr("_CoordinateAxisType","string","GeoY");
1748 
1749  at = das.get_table(dim1name);
1750  if (!at)
1751  at = das.add_table(dim1name, new AttrTable);
1752 
1753  at->append_attr("standard_name","String","projection_x_coordinate");
1754  long_name="x coordinate of projection for grid "+ gdset->getName();
1755  at->append_attr("long_name","String",long_name);
1756 
1757  // change this to meter.
1758  at->append_attr("units","string","meter");
1759  at->append_attr("_CoordinateAxisType","string","GeoX");
1760 
1761  // Add the attributes for the dummy projection variable.
1762  string cf_projection_base = "eos_cf_projection";
1763  string cf_projection = HDFCFUtil::get_CF_string(gdset->getName()) +"_"+cf_projection_base;
1764  at = das.get_table(cf_projection);
1765  if (!at)
1766  at = das.add_table(cf_projection, new AttrTable);
1767 
1768  at->append_attr("grid_mapping_name","String","sinusoidal");
1769  at->append_attr("longitude_of_central_meridian","Float64","0.0");
1770  at->append_attr("earth_radius","Float64","6371007.181");
1771 
1772  at->append_attr("_CoordinateAxisTypes","string","GeoX GeoY");
1773 
1774  // Fill in the data fields that contains the dim0name and dim1name dimensions with the grid_mapping
1775  // We only apply to >=2D data fields.
1776  HDFCFUtil::add_cf_grid_mapping_attr(das,gdset,cf_projection,dim0name,dim0size,dim1name,dim1size);
1777  }
1778 
1779 }
1780 
1781 // This function adds the 1-D horizontal coordinate variables as well as the dummy projection variable to the grid.
1782 //Note: Since we don't add these artifical CF variables to our main engineering at HDFEOS2.cc, the information
1783 // to handle DAS won't pass to DDS by the file pointer, we need to re-call the routines to check projection
1784 // and dimension. The time to retrieve these information is trivial compared with the whole translation.
1785 void HDFCFUtil::add_cf_grid_cvs(DDS & dds, HDFEOS2::GridDataset *gdset) {
1786 
1787  //1. Check the projection information, now, we only handle sinusoidal now
1788  if(GCTP_SNSOID == gdset->getProjection().getCode()) {
1789 
1790  //2. Obtain the dimension information from latitude and longitude(fieldtype =1 or fieldtype =2)
1791  string dim0name,dim1name;
1792  int32 dim0size = -1,dim1size = -1;
1793  HDFCFUtil::obtain_grid_latlon_dim_info(gdset,dim0name,dim0size,dim1name,dim1size);
1794 
1795  //3. Add the 1-D CV variables and the dummy projection variable
1796  // Note: we just need to pass the parameters that calculate 1-D cv to the data reading function,
1797  // in that way, we save the open cost of HDF-EOS2.
1798  BaseType *bt_dim0 = NULL;
1799  BaseType *bt_dim1 = NULL;
1800 
1801  HDFEOS2GeoCF1D * ar_dim0 = NULL;
1802  HDFEOS2GeoCF1D * ar_dim1 = NULL;
1803 
1804  float64 *upleft = NULL;
1805  float64 *lowright = NULL;
1806 
1807  try {
1808 
1809  bt_dim0 = new(HDFFloat64)(dim0name,gdset->getName());
1810  bt_dim1 = new(HDFFloat64)(dim1name,gdset->getName());
1811 
1812  // Obtain the upleft and lowright coordinates
1813  upleft = const_cast<float64 *>(gdset->getInfo().getUpLeft());
1814  lowright = const_cast<float64 *>(gdset->getInfo().getLowRight());
1815 
1816  // Note ar_dim0 is y, ar_dim1 is x.
1817  ar_dim0 = new HDFEOS2GeoCF1D(GCTP_SNSOID,
1818  upleft[1],lowright[1],dim0size,dim0name,bt_dim0);
1819  ar_dim0->append_dim(dim0size,dim0name);
1820 
1821  ar_dim1 = new HDFEOS2GeoCF1D(GCTP_SNSOID,
1822  upleft[0],lowright[0],dim1size,dim1name,bt_dim1);
1823  ar_dim1->append_dim(dim1size,dim1name);
1824  dds.add_var(ar_dim0);
1825  dds.add_var(ar_dim1);
1826 
1827  }
1828  catch(...) {
1829  if(bt_dim0)
1830  delete bt_dim0;
1831  if(bt_dim1)
1832  delete bt_dim1;
1833  if(ar_dim0)
1834  delete ar_dim0;
1835  if(ar_dim1)
1836  delete ar_dim1;
1837  throw InternalErr(__FILE__,__LINE__,"Unable to allocate the HDFEOS2GeoCF1D instance.");
1838  }
1839 
1840  if(bt_dim0)
1841  delete bt_dim0;
1842  if(bt_dim1)
1843  delete bt_dim1;
1844  if(ar_dim0)
1845  delete ar_dim0;
1846  if(ar_dim1)
1847  delete ar_dim1;
1848 
1849  // Also need to add the dummy projection variable.
1850  string cf_projection_base = "eos_cf_projection";
1851 
1852  // To handle multi-grid cases, we need to add the grid name.
1853  string cf_projection = HDFCFUtil::get_CF_string(gdset->getName()) +"_"+cf_projection_base;
1854 
1855  HDFEOS2GeoCFProj * dummy_proj_cf = new HDFEOS2GeoCFProj(cf_projection,gdset->getName());
1856  dds.add_var(dummy_proj_cf);
1857  if(dummy_proj_cf)
1858  delete dummy_proj_cf;
1859 
1860  }
1861 
1862 }
1863 #endif
1864 
1865  // Check OBPG attributes. Specifically, check if slope and intercept can be obtained from the file level.
1866  // If having global slope and intercept, obtain OBPG scaling, slope and intercept values.
1867 void HDFCFUtil::check_obpg_global_attrs(HDFSP::File *f, std::string & scaling,
1868  float & slope,bool &global_slope_flag,
1869  float & intercept, bool & global_intercept_flag){
1870 
1871  HDFSP::SD* spsd = f->getSD();
1872 
1873  for(vector<HDFSP::Attribute *>::const_iterator i=spsd->getAttributes().begin();i!=spsd->getAttributes().end();i++) {
1874 
1875  //We want to add two new attributes, "scale_factor" and "add_offset" to data fields if the scaling equation is linear.
1876  // OBPG products use "Slope" instead of "scale_factor", "intercept" instead of "add_offset". "Scaling" describes if the equation is linear.
1877  // Their values will be copied directly from File attributes. This addition is for OBPG L3 only.
1878  // We also add OBPG L2 support since all OBPG level 2 products(CZCS,MODISA,MODIST,OCTS,SeaWiFS) we evaluate use Slope,intercept and linear equation
1879  // for the final data. KY 2012-09-06
1880  if(f->getSPType()==OBPGL3 || f->getSPType() == OBPGL2)
1881  {
1882  if((*i)->getName()=="Scaling")
1883  {
1884  string tmpstring((*i)->getValue().begin(), (*i)->getValue().end());
1885  scaling = tmpstring;
1886  }
1887  if((*i)->getName()=="Slope" || (*i)->getName()=="slope")
1888  {
1889  global_slope_flag = true;
1890 
1891  switch((*i)->getType())
1892  {
1893 #define GET_SLOPE(TYPE, CAST) \
1894  case DFNT_##TYPE: \
1895  { \
1896  CAST tmpvalue = *(CAST*)&((*i)->getValue()[0]); \
1897  slope = (float)tmpvalue; \
1898  } \
1899  break;
1900  GET_SLOPE(INT16, int16);
1901  GET_SLOPE(INT32, int32);
1902  GET_SLOPE(FLOAT32, float);
1903  GET_SLOPE(FLOAT64, double);
1904  default:
1905  throw InternalErr(__FILE__,__LINE__,"unsupported data type.");
1906 #undef GET_SLOPE
1907  } ;
1908  //slope = *(float*)&((*i)->getValue()[0]);
1909  }
1910  if((*i)->getName()=="Intercept" || (*i)->getName()=="intercept")
1911  {
1912  global_intercept_flag = true;
1913  switch((*i)->getType())
1914  {
1915 #define GET_INTERCEPT(TYPE, CAST) \
1916  case DFNT_##TYPE: \
1917  { \
1918  CAST tmpvalue = *(CAST*)&((*i)->getValue()[0]); \
1919  intercept = (float)tmpvalue; \
1920  } \
1921  break;
1922  GET_INTERCEPT(INT16, int16);
1923  GET_INTERCEPT(INT32, int32);
1924  GET_INTERCEPT(FLOAT32, float);
1925  GET_INTERCEPT(FLOAT64, double);
1926  default:
1927  throw InternalErr(__FILE__,__LINE__,"unsupported data type.");
1928 #undef GET_INTERCEPT
1929  } ;
1930  //intercept = *(float*)&((*i)->getValue()[0]);
1931  }
1932  }
1933  }
1934 }
1935 
1936 // For some OBPG files that only provide slope and intercept at the file level,
1937 // global slope and intercept are needed to add to all fields and their names are needed to be changed to scale_factor and add_offset.
1938 // For OBPG files that provide slope and intercept at the field level, slope and intercept are needed to rename to scale_factor and add_offset.
1939 void HDFCFUtil::add_obpg_special_attrs(HDFSP::File*f,DAS &das,
1940  HDFSP::SDField *onespsds,
1941  string& scaling, float& slope,
1942  bool& global_slope_flag,
1943  float& intercept,
1944  bool & global_intercept_flag) {
1945 
1946  AttrTable *at = das.get_table(onespsds->getNewName());
1947  if (!at)
1948  at = das.add_table(onespsds->getNewName(), new AttrTable);
1949 
1950  //For OBPG L2 and L3 only.Some OBPG products put "slope" and "Intercept" etc. in the field attributes
1951  // We still need to handle the scale and offset here.
1952  bool scale_factor_flag = false;
1953  bool add_offset_flag = false;
1954  bool slope_flag = false;
1955  bool intercept_flag = false;
1956 
1957  if(f->getSPType()==OBPGL3 || f->getSPType() == OBPGL2) {// Begin OPBG CF attribute handling(Checking "slope" and "intercept")
1958  for(vector<HDFSP::Attribute *>::const_iterator i=onespsds->getAttributes().begin();
1959  i!=onespsds->getAttributes().end();i++) {
1960  if(global_slope_flag != true && ((*i)->getName()=="Slope" || (*i)->getName()=="slope"))
1961  {
1962  slope_flag = true;
1963 
1964  switch((*i)->getType())
1965  {
1966 #define GET_SLOPE(TYPE, CAST) \
1967  case DFNT_##TYPE: \
1968  { \
1969  CAST tmpvalue = *(CAST*)&((*i)->getValue()[0]); \
1970  slope = (float)tmpvalue; \
1971  } \
1972  break;
1973 
1974  GET_SLOPE(INT16, int16);
1975  GET_SLOPE(INT32, int32);
1976  GET_SLOPE(FLOAT32, float);
1977  GET_SLOPE(FLOAT64, double);
1978  default:
1979  throw InternalErr(__FILE__,__LINE__,"unsupported data type.");
1980 
1981 #undef GET_SLOPE
1982  } ;
1983  //slope = *(float*)&((*i)->getValue()[0]);
1984  }
1985  if(global_intercept_flag != true && ((*i)->getName()=="Intercept" || (*i)->getName()=="intercept"))
1986  {
1987  intercept_flag = true;
1988  switch((*i)->getType())
1989  {
1990 #define GET_INTERCEPT(TYPE, CAST) \
1991  case DFNT_##TYPE: \
1992  { \
1993  CAST tmpvalue = *(CAST*)&((*i)->getValue()[0]); \
1994  intercept = (float)tmpvalue; \
1995  } \
1996  break;
1997  GET_INTERCEPT(INT16, int16);
1998  GET_INTERCEPT(INT32, int32);
1999  GET_INTERCEPT(FLOAT32, float);
2000  GET_INTERCEPT(FLOAT64, double);
2001  default:
2002  throw InternalErr(__FILE__,__LINE__,"unsupported data type.");
2003 
2004 #undef GET_INTERCEPT
2005  } ;
2006  //intercept = *(float*)&((*i)->getValue()[0]);
2007  }
2008  }
2009  } // End of checking "slope" and "intercept"
2010 
2011  // Checking if OBPG has "scale_factor" ,"add_offset", generally checking for "long_name" attributes.
2012  for(vector<HDFSP::Attribute *>::const_iterator i=onespsds->getAttributes().begin();i!=onespsds->getAttributes().end();i++) {
2013 
2014  if((f->getSPType()==OBPGL3 || f->getSPType() == OBPGL2) && (*i)->getName()=="scale_factor")
2015  scale_factor_flag = true;
2016 
2017  if((f->getSPType()==OBPGL3 || f->getSPType() == OBPGL2) && (*i)->getName()=="add_offset")
2018  add_offset_flag = true;
2019  }
2020 
2021  // Checking if the scale and offset equation is linear, this is only for OBPG level 3.
2022  // Also checking if having the fill value and adding fill value if not having one for some OBPG data type
2023  if((f->getSPType() == OBPGL2 || f->getSPType()==OBPGL3 )&& onespsds->getFieldType()==0)
2024  {
2025 
2026  if((OBPGL3 == f->getSPType() && (scaling.find("linear")!=string::npos)) || OBPGL2 == f->getSPType() ) {
2027 
2028  if(false == scale_factor_flag && (true == slope_flag || true == global_slope_flag))
2029  {
2030  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32, 0, (void*)&(slope));
2031  at->append_attr("scale_factor", HDFCFUtil::print_type(DFNT_FLOAT32), print_rep);
2032  }
2033 
2034  if(false == add_offset_flag && (true == intercept_flag || true == global_intercept_flag))
2035  {
2036  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32, 0, (void*)&(intercept));
2037  at->append_attr("add_offset", HDFCFUtil::print_type(DFNT_FLOAT32), print_rep);
2038  }
2039  }
2040 
2041  bool has_fill_value = false;
2042  for(vector<HDFSP::Attribute *>::const_iterator i=onespsds->getAttributes().begin();i!=onespsds->getAttributes().end();i++) {
2043  if ("_FillValue" == (*i)->getNewName()){
2044  has_fill_value = true;
2045  break;
2046  }
2047  }
2048 
2049 
2050  // Add fill value to some type of OBPG data. 16-bit integer, fill value = -32767, unsigned 16-bit integer fill value = 65535
2051  // This is based on the evaluation of the example files. KY 2012-09-06
2052  if ((false == has_fill_value) &&(DFNT_INT16 == onespsds->getType())) {
2053  short fill_value = -32767;
2054  string print_rep = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)&(fill_value));
2055  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_INT16),print_rep);
2056  }
2057 
2058  if ((false == has_fill_value) &&(DFNT_UINT16 == onespsds->getType())) {
2059  unsigned short fill_value = 65535;
2060  string print_rep = HDFCFUtil::print_attr(DFNT_UINT16,0,(void*)&(fill_value));
2061  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_UINT16),print_rep);
2062  }
2063 
2064  }// Finish OBPG handling
2065 
2066 }
2067 
2068 // Handle HDF4 OTHERHDF products that follow SDS dimension scale model.
2069 // The special handling of AVHRR data is also included.
2070 void HDFCFUtil::handle_otherhdf_special_attrs(HDFSP::File*f,DAS &das) {
2071 
2072  // For some HDF4 files that follow HDF4 dimension scales, P.O. DAAC's AVHRR files.
2073  // The "otherHDF" category can almost make AVHRR files work, except
2074  // that AVHRR uses the attribute name "unit" instead of "units" for latitude and longitude,
2075  // I have to correct the name to follow CF conventions(using "units"). I won't check
2076  // the latitude and longitude values since latitude and longitude values for some files(LISO files)
2077  // are not in the standard range(0-360 for lon and 0-180 for lat). KY 2011-3-3
2078  const vector<HDFSP::SDField *>& spsds = f->getSD()->getFields();
2079  vector<HDFSP::SDField *>::const_iterator it_g;
2080 
2081  if(f->getSPType() == OTHERHDF) {
2082 
2083  bool latflag = false;
2084  bool latunitsflag = false; //CF conventions use "units" instead of "unit"
2085  bool lonflag = false;
2086  bool lonunitsflag = false; // CF conventions use "units" instead of "unit"
2087  int llcheckoverflag = 0;
2088 
2089  // Here I try to condense the code within just two for loops
2090  // The outer loop: Loop through all SDS objects
2091  // The inner loop: Loop through all attributes of this SDS
2092  // Inside two inner loops(since "units" are common for other fields),
2093  // inner loop 1: when an attribute's long name value is L(l)atitude,mark it.
2094  // inner loop 2: when an attribute's name is units, mark it.
2095  // Outside the inner loop, when latflag is true, and latunitsflag is false,
2096  // adding a new attribute called units and the value should be "degrees_north".
2097  // doing the same thing for longitude.
2098 
2099  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2100 
2101  // Ignore ALL coordinate variables if this is "OTHERHDF" case and some dimensions
2102  // don't have dimension scale data.
2103  if ( true == f->Has_Dim_NoScale_Field() && ((*it_g)->getFieldType() !=0) && ((*it_g)->IsDimScale() == false))
2104  continue;
2105 
2106  // Ignore the empty(no data) dimension variable.
2107  if (OTHERHDF == f->getSPType() && true == (*it_g)->IsDimNoScale())
2108  continue;
2109 
2110  AttrTable *at = das.get_table((*it_g)->getNewName());
2111  if (!at)
2112  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2113 
2114  for(vector<HDFSP::Attribute *>::const_iterator i=(*it_g)->getAttributes().begin();i!=(*it_g)->getAttributes().end();i++) {
2115  if((*i)->getType()==DFNT_UCHAR || (*i)->getType() == DFNT_CHAR){
2116 
2117  if((*i)->getName() == "long_name") {
2118  string tempstring2((*i)->getValue().begin(),(*i)->getValue().end());
2119  string tempfinalstr= string(tempstring2.c_str());// This may remove some garbage characters
2120  if(tempfinalstr=="latitude" || tempfinalstr == "Latitude") // Find long_name latitude
2121  latflag = true;
2122  if(tempfinalstr=="longitude" || tempfinalstr == "Longitude") // Find long_name latitude
2123  lonflag = true;
2124  }
2125  }
2126  }
2127 
2128  if(latflag) {
2129  for(vector<HDFSP::Attribute *>::const_iterator i=(*it_g)->getAttributes().begin();i!=(*it_g)->getAttributes().end();i++) {
2130  if((*i)->getName() == "units")
2131  latunitsflag = true;
2132  }
2133  }
2134 
2135  if(lonflag) {
2136  for(vector<HDFSP::Attribute *>::const_iterator i=(*it_g)->getAttributes().begin();i!=(*it_g)->getAttributes().end();i++) {
2137  if((*i)->getName() == "units")
2138  lonunitsflag = true;
2139  }
2140  }
2141  if(latflag && !latunitsflag){ // No "units" for latitude, add "units"
2142  at->append_attr("units","String","degrees_north");
2143  latflag = false;
2144  latunitsflag = false;
2145  llcheckoverflag++;
2146  }
2147 
2148  if(lonflag && !lonunitsflag){ // No "units" for latitude, add "units"
2149  at->append_attr("units","String","degrees_east");
2150  lonflag = false;
2151  lonunitsflag = false;
2152  llcheckoverflag++;
2153  }
2154  if(llcheckoverflag ==2) break;
2155 
2156  }
2157 
2158  }
2159 
2160 }
2161 
2162 // Add missing CF attributes for non-CV varibles
2163 void
2164 HDFCFUtil::add_missing_cf_attrs(HDFSP::File*f,DAS &das) {
2165 
2166 
2167  const vector<HDFSP::SDField *>& spsds = f->getSD()->getFields();
2168  vector<HDFSP::SDField *>::const_iterator it_g;
2169 
2170 
2171  // TRMM level 3 grid
2172  if(TRMML3A_V6== f->getSPType() || TRMML3C_V6==f->getSPType() || TRMML3S_V7 == f->getSPType() || TRMML3M_V7 == f->getSPType()) {
2173 
2174  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2175  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_FLOAT32) {
2176 
2177  AttrTable *at = das.get_table((*it_g)->getNewName());
2178  if (!at)
2179  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2180  string print_rep = "-9999.9";
2181  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_FLOAT32),print_rep);
2182 
2183  }
2184  }
2185 
2186  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2187  if((*it_g)->getFieldType() == 0 && (((*it_g)->getType()==DFNT_INT32) || ((*it_g)->getType()==DFNT_INT16))) {
2188 
2189  AttrTable *at = das.get_table((*it_g)->getNewName());
2190  if (!at)
2191  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2192  string print_rep = "-9999";
2193  if((*it_g)->getType()==DFNT_INT32)
2194  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_INT32),print_rep);
2195  else
2196  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_INT16),print_rep);
2197 
2198  }
2199  }
2200 
2201  // nlayer for TRMM single grid version 7, the units should be "km"
2202  if(TRMML3S_V7 == f->getSPType()) {
2203  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2204  if((*it_g)->getFieldType() == 6 && (*it_g)->getNewName()=="nlayer") {
2205 
2206  AttrTable *at = das.get_table((*it_g)->getNewName());
2207  if (!at)
2208  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2209  at->append_attr("units","String","km");
2210 
2211  }
2212  else if((*it_g)->getFieldType() == 4) {
2213 
2214  if ((*it_g)->getNewName()=="nh3" ||
2215  (*it_g)->getNewName()=="ncat3" ||
2216  (*it_g)->getNewName()=="nthrshZO" ||
2217  (*it_g)->getNewName()=="nthrshHB" ||
2218  (*it_g)->getNewName()=="nthrshSRT")
2219  {
2220 
2221  string references =
2222  "http://pps.gsfc.nasa.gov/Documents/filespec.TRMM.V7.pdf";
2223  string comment;
2224 
2225  AttrTable *at = das.get_table((*it_g)->getNewName());
2226  if (!at)
2227  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2228 
2229  if((*it_g)->getNewName()=="nh3") {
2230  comment="Index number to represent the fixed heights above the earth ellipsoid,";
2231  comment= comment + " at 2, 4, 6 km plus one for path-average.";
2232  }
2233 
2234  else if((*it_g)->getNewName()=="ncat3") {
2235  comment="Index number to represent catgories for probability distribution functions.";
2236  comment=comment + "Check more information from the references.";
2237  }
2238 
2239  else if((*it_g)->getNewName()=="nthrshZO")
2240  comment="Q-thresholds for Zero order used for probability distribution functions.";
2241 
2242  else if((*it_g)->getNewName()=="nthrshHB")
2243  comment="Q-thresholds for HB used for probability distribution functions.";
2244 
2245  else if((*it_g)->getNewName()=="nthrshSRT")
2246  comment="Q-thresholds for SRT used for probability distribution functions.";
2247 
2248  at->append_attr("comment","String",comment);
2249  at->append_attr("references","String",references);
2250 
2251  }
2252 
2253  }
2254 
2255  }
2256 
2257 
2258  // 3A26 use special values such as -666, -777,-999 in their fields.
2259  // Although the document doesn't provide range for some fields, the meaning of those fields should be greater than 0.
2260  // So add valid_min = 0 and fill_value = -999 .
2261  string base_filename;
2262  size_t last_slash_pos = f->getPath().find_last_of("/");
2263  if(last_slash_pos != string::npos)
2264  base_filename = f->getPath().substr(last_slash_pos+1);
2265  if(""==base_filename)
2266  base_filename = f->getPath();
2267  bool t3a26_flag = ((base_filename.find("3A26")!=string::npos)?true:false);
2268 
2269  if(true == t3a26_flag) {
2270 
2271  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2272 
2273  if((*it_g)->getFieldType() == 0 && ((*it_g)->getType()==DFNT_FLOAT32)) {
2274  AttrTable *at = das.get_table((*it_g)->getNewName());
2275  if (!at)
2276  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2277  at->del_attr("_FillValue");
2278  at->append_attr("_FillValue","Float32","-999");
2279  at->append_attr("valid_min","Float32","0");
2280 
2281  }
2282  }
2283  }
2284 
2285  }
2286 
2287  // nlayer for TRMM single grid version 7, the units should be "km"
2288  if(TRMML3M_V7 == f->getSPType()) {
2289  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2290 
2291  if((*it_g)->getFieldType() == 4 ) {
2292 
2293  string references ="http://pps.gsfc.nasa.gov/Documents/filespec.TRMM.V7.pdf";
2294  if ((*it_g)->getNewName()=="nh1") {
2295 
2296  AttrTable *at = das.get_table((*it_g)->getNewName());
2297  if (!at)
2298  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2299 
2300  string comment="Number of fixed heights above the earth ellipsoid,";
2301  comment= comment + " at 2, 4, 6, 10, and 15 km plus one for path-average.";
2302 
2303  at->append_attr("comment","String",comment);
2304  at->append_attr("references","String",references);
2305 
2306  }
2307  if ((*it_g)->getNewName()=="nh3") {
2308 
2309  AttrTable *at = das.get_table((*it_g)->getNewName());
2310  if (!at)
2311  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2312 
2313  string comment="Number of fixed heights above the earth ellipsoid,";
2314  comment= comment + " at 2, 4, 6 km plus one for path-average.";
2315 
2316  at->append_attr("comment","String",comment);
2317  at->append_attr("references","String",references);
2318 
2319  }
2320 
2321  if ((*it_g)->getNewName()=="nang") {
2322 
2323  AttrTable *at = das.get_table((*it_g)->getNewName());
2324  if (!at)
2325  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2326 
2327  string comment="Number of fixed incidence angles, at 0, 5, 10 and 15 degree and all angles.";
2328  references = "http://pps.gsfc.nasa.gov/Documents/ICSVol4.pdf";
2329 
2330  at->append_attr("comment","String",comment);
2331  at->append_attr("references","String",references);
2332 
2333  }
2334 
2335  if ((*it_g)->getNewName()=="ncat2") {
2336 
2337  AttrTable *at = das.get_table((*it_g)->getNewName());
2338  if (!at)
2339  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2340 
2341  string comment="Second number of categories for histograms (30). ";
2342  comment=comment + "Check more information from the references.";
2343 
2344  at->append_attr("comment","String",comment);
2345  at->append_attr("references","String",references);
2346 
2347  }
2348 
2349  }
2350  }
2351 
2352  }
2353 
2354  }
2355 
2356  // TRMM level 2 swath
2357  else if(TRMML2_V7== f->getSPType()) {
2358 
2359  string base_filename;
2360  size_t last_slash_pos = f->getPath().find_last_of("/");
2361  if(last_slash_pos != string::npos)
2362  base_filename = f->getPath().substr(last_slash_pos+1);
2363  if(""==base_filename)
2364  base_filename = f->getPath();
2365  bool t2b31_flag = ((base_filename.find("2B31")!=string::npos)?true:false);
2366  bool t2a21_flag = ((base_filename.find("2A21")!=string::npos)?true:false);
2367  bool t2a12_flag = ((base_filename.find("2A12")!=string::npos)?true:false);
2368  // 2A23 is temporarily not supported perhaps due to special fill values
2369  //bool t2a23_flag = ((base_filename.find("2A23")!=string::npos)?true:false);
2370  bool t2a25_flag = ((base_filename.find("2A25")!=string::npos)?true:false);
2371  bool t1c21_flag = ((base_filename.find("1C21")!=string::npos)?true:false);
2372  bool t1b21_flag = ((base_filename.find("1B21")!=string::npos)?true:false);
2373  bool t1b11_flag = ((base_filename.find("1B11")!=string::npos)?true:false);
2374  bool t1b01_flag = ((base_filename.find("1B01")!=string::npos)?true:false);
2375 
2376  // Handle scale and offset
2377 
2378  // group 1: 2B31,2A12,2A21
2379  if(t2b31_flag || t2a12_flag || t2a21_flag) {
2380 
2381  // special for 2B31
2382  if(t2b31_flag) {
2383 
2384  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2385 
2386 
2387  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_INT16) {
2388 
2389  AttrTable *at = das.get_table((*it_g)->getNewName());
2390  if (!at)
2391  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2392 
2393  AttrTable::Attr_iter it = at->attr_begin();
2394  while (it!=at->attr_end()) {
2395  if(at->get_name(it)=="scale_factor")
2396  {
2397  // Scale type and value
2398  string scale_factor_value="";
2399  string scale_factor_type;
2400 
2401  scale_factor_value = (*at->get_attr_vector(it)->begin());
2402  scale_factor_type = at->get_type(it);
2403 
2404  if(scale_factor_type == "Float64") {
2405  double new_scale = 1.0/strtod(scale_factor_value.c_str(),NULL);
2406  at->del_attr("scale_factor");
2407  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT64,0,(void*)(&new_scale));
2408  at->append_attr("scale_factor", scale_factor_type,print_rep);
2409 
2410  }
2411 
2412  if(scale_factor_type == "Float32") {
2413  float new_scale = 1.0/strtof(scale_factor_value.c_str(),NULL);
2414  at->del_attr("scale_factor");
2415  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&new_scale));
2416  at->append_attr("scale_factor", scale_factor_type,print_rep);
2417 
2418  }
2419 
2420  break;
2421  }
2422  ++it;
2423  }
2424  }
2425  }
2426  }
2427 
2428  // Special for 2A12
2429  if(t2a12_flag==true) {
2430 
2431  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2432 
2433  if((*it_g)->getFieldType() == 6 && (*it_g)->getNewName()=="nlayer") {
2434 
2435  AttrTable *at = das.get_table((*it_g)->getNewName());
2436  if (!at)
2437  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2438  at->append_attr("units","String","km");
2439 
2440  }
2441 
2442  // signed char maps to int32, so use int32 for the fillvalue.
2443  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_INT8) {
2444 
2445  AttrTable *at = das.get_table((*it_g)->getNewName());
2446  if (!at)
2447  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2448  at->append_attr("_FillValue","Int32","-99");
2449 
2450  }
2451 
2452  }
2453  }
2454 
2455 
2456  // for all 2A12,2A21 and 2B31
2457  // Add fillvalues for float32 and int32.
2458  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2459  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_FLOAT32) {
2460 
2461  AttrTable *at = das.get_table((*it_g)->getNewName());
2462  if (!at)
2463  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2464  string print_rep = "-9999.9";
2465  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_FLOAT32),print_rep);
2466 
2467  }
2468  }
2469 
2470  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2471 
2472 
2473  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_INT16) {
2474 
2475  AttrTable *at = das.get_table((*it_g)->getNewName());
2476  if (!at)
2477  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2478 
2479  string print_rep = "-9999";
2480  at->append_attr("_FillValue",HDFCFUtil::print_type(DFNT_INT32),print_rep);
2481 
2482  }
2483  }
2484 
2485  }
2486 
2487  // group 2: 2A21 and 2A25.
2488  else if(t2a21_flag == true || t2a25_flag == true) {
2489 
2490  // 2A25: handle reflectivity and rain rate scales
2491  if(t2a25_flag == true) {
2492 
2493  unsigned char handle_scale = 0;
2494  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2495 
2496  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_INT16) {
2497  bool has_dBZ = false;
2498  bool has_rainrate = false;
2499  bool has_scale = false;
2500  string scale_factor_value;
2501  string scale_factor_type;
2502 
2503  AttrTable *at = das.get_table((*it_g)->getNewName());
2504  if (!at)
2505  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2506  AttrTable::Attr_iter it = at->attr_begin();
2507  while (it!=at->attr_end()) {
2508  if(at->get_name(it)=="units"){
2509  string units_value = *at->get_attr_vector(it)->begin();
2510  if("dBZ" == units_value) {
2511  has_dBZ = true;
2512  }
2513 
2514  else if("mm/hr" == units_value){
2515  has_rainrate = true;
2516  }
2517  }
2518  if(at->get_name(it)=="scale_factor")
2519  {
2520  scale_factor_value = *at->get_attr_vector(it)->begin();
2521  scale_factor_type = at->get_type(it);
2522  has_scale = true;
2523  }
2524  ++it;
2525 
2526  }
2527 
2528  if((true == has_rainrate || true == has_dBZ) && true == has_scale) {
2529 
2530  handle_scale++;
2531  short valid_min = 0;
2532  short valid_max = 0;
2533 
2534  // Here just use 32-bit floating-point for the scale_factor, should be okay.
2535  if(true == has_rainrate)
2536  valid_max = (short)(300*strtof(scale_factor_value.c_str(),NULL));
2537  else if(true == has_dBZ)
2538  valid_max = (short)(80*strtof(scale_factor_value.c_str(),NULL));
2539 
2540  string print_rep1 = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)(&valid_min));
2541  at->append_attr("valid_min","Int16",print_rep1);
2542  print_rep1 = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)(&valid_max));
2543  at->append_attr("valid_max","Int16",print_rep1);
2544 
2545  at->del_attr("scale_factor");
2546  if(scale_factor_type == "Float64") {
2547  double new_scale = 1.0/strtod(scale_factor_value.c_str(),NULL);
2548  string print_rep2 = HDFCFUtil::print_attr(DFNT_FLOAT64,0,(void*)(&new_scale));
2549  at->append_attr("scale_factor", scale_factor_type,print_rep2);
2550 
2551  }
2552 
2553  if(scale_factor_type == "Float32") {
2554  float new_scale = 1.0/strtof(scale_factor_value.c_str(),NULL);
2555  string print_rep3 = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&new_scale));
2556  at->append_attr("scale_factor", scale_factor_type,print_rep3);
2557 
2558  }
2559 
2560 
2561  }
2562 
2563  if(2 == handle_scale)
2564  break;
2565 
2566  }
2567  }
2568  }
2569  }
2570 
2571  // 1B21,1C21 and 1B11
2572  else if(t1b21_flag || t1c21_flag || t1b11_flag) {
2573 
2574  // 1B21,1C21 scale_factor to CF and valid_range for dBm and dBZ.
2575  if(t1b21_flag || t1c21_flag) {
2576 
2577  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2578 
2579  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_INT16) {
2580 
2581  bool has_dBm = false;
2582  bool has_dBZ = false;
2583 
2584  AttrTable *at = das.get_table((*it_g)->getNewName());
2585  if (!at)
2586  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2587  AttrTable::Attr_iter it = at->attr_begin();
2588 
2589  while (it!=at->attr_end()) {
2590  if(at->get_name(it)=="units"){
2591 
2592  string units_value = *at->get_attr_vector(it)->begin();
2593  if("dBm" == units_value) {
2594  has_dBm = true;
2595  break;
2596  }
2597 
2598  else if("dBZ" == units_value){
2599  has_dBZ = true;
2600  break;
2601  }
2602  }
2603  ++it;
2604  }
2605 
2606  if(has_dBm == true || has_dBZ == true) {
2607  it = at->attr_begin();
2608  while (it!=at->attr_end()) {
2609  if(at->get_name(it)=="scale_factor")
2610  {
2611 
2612  string scale_value = *at->get_attr_vector(it)->begin();
2613 
2614  if(true == has_dBm) {
2615  short valid_min = (short)(-120 *strtof(scale_value.c_str(),NULL));
2616  short valid_max = (short)(-20 *strtof(scale_value.c_str(),NULL));
2617  string print_rep = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)(&valid_min));
2618  at->append_attr("valid_min","Int16",print_rep);
2619  print_rep = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)(&valid_max));
2620  at->append_attr("valid_max","Int16",print_rep);
2621  break;
2622 
2623  }
2624 
2625  else if(true == has_dBZ){
2626  short valid_min = (short)(-20 *strtof(scale_value.c_str(),NULL));
2627  short valid_max = (short)(80 *strtof(scale_value.c_str(),NULL));
2628  string print_rep = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)(&valid_min));
2629  at->append_attr("valid_min","Int16",print_rep);
2630  print_rep = HDFCFUtil::print_attr(DFNT_INT16,0,(void*)(&valid_max));
2631  at->append_attr("valid_max","Int16",print_rep);
2632  break;
2633 
2634  }
2635  }
2636  ++it;
2637  }
2638 
2639  }
2640  }
2641  }
2642  }
2643 
2644  // For all 1B21,1C21 and 1B11 int16-bit products,change scale to follow CF
2645  // I find that one 1B21 variable binStormHeight has fillvalue -9999,
2646  // so add _FillValue -9999 for int16-bit variables.
2647  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2648 
2649  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_INT16) {
2650 
2651  AttrTable *at = das.get_table((*it_g)->getNewName());
2652  if (!at)
2653  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2654  AttrTable::Attr_iter it = at->attr_begin();
2655 
2656 
2657  while (it!=at->attr_end()) {
2658 
2659  if(at->get_name(it)=="scale_factor")
2660  {
2661  // Scale type and value
2662  string scale_factor_value="";
2663  string scale_factor_type;
2664 
2665  scale_factor_value = (*at->get_attr_vector(it)->begin());
2666  scale_factor_type = at->get_type(it);
2667 
2668  if(scale_factor_type == "Float64") {
2669  double new_scale = 1.0/strtod(scale_factor_value.c_str(),NULL);
2670  at->del_attr("scale_factor");
2671  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT64,0,(void*)(&new_scale));
2672  at->append_attr("scale_factor", scale_factor_type,print_rep);
2673 
2674  }
2675 
2676  if(scale_factor_type == "Float32") {
2677  float new_scale = 1.0/strtof(scale_factor_value.c_str(),NULL);
2678  at->del_attr("scale_factor");
2679  string print_rep = HDFCFUtil::print_attr(DFNT_FLOAT32,0,(void*)(&new_scale));
2680  at->append_attr("scale_factor", scale_factor_type,print_rep);
2681 
2682  }
2683 
2684  break;
2685 
2686  }
2687  ++it;
2688 
2689  }
2690 
2691  at->append_attr("_FillValue","Int16","-9999");
2692 
2693  }
2694  }
2695  }
2696 
2697  // For 1B01 product, just add the fillvalue.
2698  else if(t1b01_flag == true) {
2699  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2700 
2701  if((*it_g)->getFieldType() == 0 && (*it_g)->getType()==DFNT_FLOAT32) {
2702  AttrTable *at = das.get_table((*it_g)->getNewName());
2703  if (!at)
2704  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2705 
2706  at->append_attr("_FillValue","Float32","-9999.9");
2707 
2708  }
2709  }
2710 
2711  }
2712 
2713  AttrTable *at = das.get_table("HDF_GLOBAL");
2714  if (!at)
2715  at = das.add_table("HDF_GLOBAL", new AttrTable);
2716  string references ="http://pps.gsfc.nasa.gov/Documents/filespec.TRMM.V7.pdf";
2717  string comment="The HDF4 OPeNDAP handler adds _FillValue, valid_min and valid_max for some TRMM level 1 and level 2 products.";
2718  comment= comment + " It also changes scale_factor to follow CF conventions. ";
2719 
2720  at->append_attr("comment","String",comment);
2721  at->append_attr("references","String",references);
2722 
2723  }
2724 
2725 }
2726 
2727 
2728 
2729 //
2730 // Many CERES products compose of multiple groups
2731 // There are many fields in CERES data(a few hundred) and the full name(with the additional path)
2732 // is very long. It causes Java clients choken since Java clients append names in the URL
2733 // To improve the performance and to make Java clients access the data, we will ignore
2734 // the path of these fields. Users can turn off this feature by commenting out the line: H4.EnableCERESMERRAShortName=true
2735 // or can set the H4.EnableCERESMERRAShortName=false
2736 // We still preserve the full path as an attribute in case users need to check them.
2737 // Kent 2012-6-29
2738 
2739 void HDFCFUtil::handle_merra_ceres_attrs_with_bes_keys(HDFSP::File *f, DAS &das,const string& filename) {
2740 
2741  string base_filename = filename.substr(filename.find_last_of("/")+1);
2742 
2743 #if 0
2744  string check_ceres_merra_short_name_key="H4.EnableCERESMERRAShortName";
2745  bool turn_on_ceres_merra_short_name_key= false;
2746 
2747  turn_on_ceres_merra_short_name_key = HDFCFUtil::check_beskeys(check_ceres_merra_short_name_key);
2748 #endif
2749 
2750  bool merra_is_eos2 = false;
2751  if(0== (base_filename.compare(0,5,"MERRA"))) {
2752 
2753  for (vector < HDFSP::Attribute * >::const_iterator i =
2754  f->getSD()->getAttributes ().begin ();
2755  i != f->getSD()->getAttributes ().end (); ++i) {
2756 
2757  // CHeck if this MERRA file is an HDF-EOS2 or not.
2758  if(((*i)->getName().compare(0, 14, "StructMetadata" )== 0) ||
2759  ((*i)->getName().compare(0, 14, "structmetadata" )== 0)) {
2760  merra_is_eos2 = true;
2761  break;
2762  }
2763 
2764  }
2765  }
2766 
2767  if (true == HDF4RequestHandler::get_enable_ceres_merra_short_name() && (CER_ES4 == f->getSPType() || CER_SRB == f->getSPType()
2768  || CER_CDAY == f->getSPType() || CER_CGEO == f->getSPType()
2769  || CER_SYN == f->getSPType() || CER_ZAVG == f->getSPType()
2770  || CER_AVG == f->getSPType() || (true == merra_is_eos2))) {
2771 
2772  const vector<HDFSP::SDField *>& spsds = f->getSD()->getFields();
2773  vector<HDFSP::SDField *>::const_iterator it_g;
2774  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
2775 
2776  AttrTable *at = das.get_table((*it_g)->getNewName());
2777  if (!at)
2778  at = das.add_table((*it_g)->getNewName(), new AttrTable);
2779 
2780  at->append_attr("fullpath","String",(*it_g)->getSpecFullPath());
2781 
2782  }
2783 
2784  }
2785 
2786 }
2787 
2788 
2789 // Handle the attributes when the BES key EnableVdataDescAttr is enabled..
2790 void HDFCFUtil::handle_vdata_attrs_with_desc_key(HDFSP::File*f,libdap::DAS &das) {
2791 
2792  // Check the EnableVdataDescAttr key. If this key is turned on, the handler-added attribute VDdescname and
2793  // the attributes of vdata and vdata fields will be outputed to DAS. Otherwise, these attributes will
2794  // not outputed to DAS. The key will be turned off by default to shorten the DAP output. KY 2012-09-18
2795 
2796 #if 0
2797  string check_vdata_desc_key="H4.EnableVdataDescAttr";
2798  bool turn_on_vdata_desc_key= false;
2799 
2800  turn_on_vdata_desc_key = HDFCFUtil::check_beskeys(check_vdata_desc_key);
2801 #endif
2802 
2803  string VDdescname = "hdf4_vd_desc";
2804  string VDdescvalue = "This is an HDF4 Vdata.";
2805  string VDfieldprefix = "Vdata_field_";
2806  string VDattrprefix = "Vdata_attr_";
2807  string VDfieldattrprefix ="Vdata_field_attr_";
2808 
2809  // To speed up the performance for handling CERES data, we turn off some CERES vdata fields, this should be resumed in the future version with BESKeys.
2810 #if 0
2811  string check_ceres_vdata_key="H4.EnableCERESVdata";
2812  bool turn_on_ceres_vdata_key= false;
2813  turn_on_ceres_vdata_key = HDFCFUtil::check_beskeys(check_ceres_vdata_key);
2814 #endif
2815 
2816  bool output_vdata_flag = true;
2817  if (false == HDF4RequestHandler::get_enable_ceres_vdata() &&
2818  (CER_AVG == f->getSPType() ||
2819  CER_ES4 == f->getSPType() ||
2820  CER_SRB == f->getSPType() ||
2821  CER_ZAVG == f->getSPType()))
2822  output_vdata_flag = false;
2823 
2824  if (true == output_vdata_flag) {
2825 
2826  for(vector<HDFSP::VDATA *>::const_iterator i=f->getVDATAs().begin(); i!=f->getVDATAs().end();i++) {
2827 
2828  AttrTable *at = das.get_table((*i)->getNewName());
2829  if(!at)
2830  at = das.add_table((*i)->getNewName(),new AttrTable);
2831 
2832  if (true == HDF4RequestHandler::get_enable_vdata_desc_attr()) {
2833  // Add special vdata attributes
2834  bool emptyvddasflag = true;
2835  if(!((*i)->getAttributes().empty())) emptyvddasflag = false;
2836  if((*i)->getTreatAsAttrFlag())
2837  emptyvddasflag = false;
2838  else {
2839  for(vector<HDFSP::VDField *>::const_iterator j=(*i)->getFields().begin();j!=(*i)->getFields().end();j++) {
2840  if(!((*j)->getAttributes().empty())) {
2841  emptyvddasflag = false;
2842  break;
2843  }
2844  }
2845  }
2846 
2847  if(emptyvddasflag)
2848  continue;
2849  at->append_attr(VDdescname, "String" , VDdescvalue);
2850 
2851  for(vector<HDFSP::Attribute *>::const_iterator it_va = (*i)->getAttributes().begin();it_va!=(*i)->getAttributes().end();it_va++) {
2852 
2853  if((*it_va)->getType()==DFNT_UCHAR || (*it_va)->getType() == DFNT_CHAR){
2854 
2855  string tempstring2((*it_va)->getValue().begin(),(*it_va)->getValue().end());
2856  string tempfinalstr= string(tempstring2.c_str());
2857  at->append_attr(VDattrprefix+(*it_va)->getNewName(), "String" , HDFCFUtil::escattr(tempfinalstr));
2858  }
2859  else {
2860  for (int loc=0; loc < (*it_va)->getCount() ; loc++) {
2861  string print_rep = HDFCFUtil::print_attr((*it_va)->getType(), loc, (void*) &((*it_va)->getValue()[0]));
2862  at->append_attr(VDattrprefix+(*it_va)->getNewName(), HDFCFUtil::print_type((*it_va)->getType()), print_rep);
2863  }
2864  }
2865 
2866  }
2867  }
2868 
2869  if(false == ((*i)->getTreatAsAttrFlag())){
2870 
2871  if (true == HDF4RequestHandler::get_enable_vdata_desc_attr()) {
2872 
2873  //NOTE: for vdata field, we assume that no special characters are found. We need to escape the special characters when the data type is char.
2874  // We need to create a DAS container for each field so that the attributes can be put inside.
2875  for(vector<HDFSP::VDField *>::const_iterator j=(*i)->getFields().begin();j!=(*i)->getFields().end();j++) {
2876 
2877  // This vdata field will NOT be treated as attributes, only save the field attribute as the attribute
2878  // First check if the field has attributes, if it doesn't have attributes, no need to create a container.
2879 
2880  if((*j)->getAttributes().size() !=0) {
2881 
2882  AttrTable *at_v = das.get_table((*j)->getNewName());
2883  if(!at_v)
2884  at_v = das.add_table((*j)->getNewName(),new AttrTable);
2885 
2886  for(vector<HDFSP::Attribute *>::const_iterator it_va = (*j)->getAttributes().begin();it_va!=(*j)->getAttributes().end();it_va++) {
2887 
2888  if((*it_va)->getType()==DFNT_UCHAR || (*it_va)->getType() == DFNT_CHAR){
2889 
2890  string tempstring2((*it_va)->getValue().begin(),(*it_va)->getValue().end());
2891  string tempfinalstr= string(tempstring2.c_str());
2892  at_v->append_attr((*it_va)->getNewName(), "String" , HDFCFUtil::escattr(tempfinalstr));
2893  }
2894  else {
2895  for (int loc=0; loc < (*it_va)->getCount() ; loc++) {
2896  string print_rep = HDFCFUtil::print_attr((*it_va)->getType(), loc, (void*) &((*it_va)->getValue()[0]));
2897  at_v->append_attr((*it_va)->getNewName(), HDFCFUtil::print_type((*it_va)->getType()), print_rep);
2898  }
2899  }
2900 
2901  }
2902  }
2903  }
2904  }
2905 
2906  }
2907 
2908  else {
2909 
2910  for(vector<HDFSP::VDField *>::const_iterator j=(*i)->getFields().begin();j!=(*i)->getFields().end();j++) {
2911 
2912  if((*j)->getFieldOrder() == 1) {
2913  if((*j)->getType()==DFNT_UCHAR || (*j)->getType() == DFNT_CHAR){
2914  string tempfinalstr;
2915  tempfinalstr.resize((*j)->getValue().size());
2916  copy((*j)->getValue().begin(),(*j)->getValue().end(),tempfinalstr.begin());
2917  at->append_attr(VDfieldprefix+(*j)->getNewName(), "String" , HDFCFUtil::escattr(tempfinalstr));
2918  }
2919  else {
2920  for ( int loc=0; loc < (*j)->getNumRec(); loc++) {
2921  string print_rep = HDFCFUtil::print_attr((*j)->getType(), loc, (void*) &((*j)->getValue()[0]));
2922  at->append_attr(VDfieldprefix+(*j)->getNewName(), HDFCFUtil::print_type((*j)->getType()), print_rep);
2923  }
2924  }
2925  }
2926  else {//When field order is greater than 1,we want to print each record in group with single quote,'0 1 2','3 4 5', etc.
2927 
2928  if((*j)->getValue().size() != (unsigned int)(DFKNTsize((*j)->getType())*((*j)->getFieldOrder())*((*j)->getNumRec()))){
2929  throw InternalErr(__FILE__,__LINE__,"the vdata field size doesn't match the vector value");
2930  }
2931 
2932  if((*j)->getNumRec()==1){
2933  if((*j)->getType()==DFNT_UCHAR || (*j)->getType() == DFNT_CHAR){
2934  string tempstring2((*j)->getValue().begin(),(*j)->getValue().end());
2935  string tempfinalstr= string(tempstring2.c_str());
2936  at->append_attr(VDfieldprefix+(*j)->getNewName(),"String",HDFCFUtil::escattr(tempfinalstr));
2937  }
2938  else {
2939  for (int loc=0; loc < (*j)->getFieldOrder(); loc++) {
2940  string print_rep = HDFCFUtil::print_attr((*j)->getType(), loc, (void*) &((*j)->getValue()[0]));
2941  at->append_attr(VDfieldprefix+(*j)->getNewName(), HDFCFUtil::print_type((*j)->getType()), print_rep);
2942  }
2943  }
2944 
2945  }
2946 
2947  else {
2948  if((*j)->getType()==DFNT_UCHAR || (*j)->getType() == DFNT_CHAR){
2949 
2950  for(int tempcount = 0; tempcount < (*j)->getNumRec()*DFKNTsize((*j)->getType());tempcount ++) {
2951  vector<char>::const_iterator tempit;
2952  tempit = (*j)->getValue().begin()+tempcount*((*j)->getFieldOrder());
2953  string tempstring2(tempit,tempit+(*j)->getFieldOrder());
2954  string tempfinalstr= string(tempstring2.c_str());
2955  string tempoutstring = "'"+tempfinalstr+"'";
2956  at->append_attr(VDfieldprefix+(*j)->getNewName(),"String",HDFCFUtil::escattr(tempoutstring));
2957  }
2958  }
2959 
2960  else {
2961  for(int tempcount = 0; tempcount < (*j)->getNumRec();tempcount ++) {
2962  at->append_attr(VDfieldprefix+(*j)->getNewName(),HDFCFUtil::print_type((*j)->getType()),"'");
2963  for (int loc=0; loc < (*j)->getFieldOrder(); loc++) {
2964  string print_rep = HDFCFUtil::print_attr((*j)->getType(), loc, (void*) &((*j)->getValue()[tempcount*((*j)->getFieldOrder())]));
2965  at->append_attr(VDfieldprefix+(*j)->getNewName(), HDFCFUtil::print_type((*j)->getType()), print_rep);
2966  }
2967  at->append_attr(VDfieldprefix+(*j)->getNewName(),HDFCFUtil::print_type((*j)->getType()),"'");
2968  }
2969  }
2970  }
2971  }
2972 
2973 
2974  if (true == HDF4RequestHandler::get_enable_vdata_desc_attr()) {
2975  for(vector<HDFSP::Attribute *>::const_iterator it_va = (*j)->getAttributes().begin();it_va!=(*j)->getAttributes().end();it_va++) {
2976 
2977  if((*it_va)->getType()==DFNT_UCHAR || (*it_va)->getType() == DFNT_CHAR){
2978 
2979  string tempstring2((*it_va)->getValue().begin(),(*it_va)->getValue().end());
2980  string tempfinalstr= string(tempstring2.c_str());
2981  at->append_attr(VDfieldattrprefix+(*it_va)->getNewName(), "String" , HDFCFUtil::escattr(tempfinalstr));
2982  }
2983  else {
2984  for (int loc=0; loc < (*it_va)->getCount() ; loc++) {
2985  string print_rep = HDFCFUtil::print_attr((*it_va)->getType(), loc, (void*) &((*it_va)->getValue()[0]));
2986  at->append_attr(VDfieldattrprefix+(*it_va)->getNewName(), HDFCFUtil::print_type((*it_va)->getType()), print_rep);
2987  }
2988  }
2989  }
2990  }
2991  }
2992  }
2993 
2994  }
2995  }
2996 
2997 }
2998 
2999 void HDFCFUtil::map_eos2_objects_attrs(libdap::DAS &das,const string &filename) {
3000 
3001  intn status_n =-1;
3002  int32 status_32 = -1;
3003  int32 file_id = -1;
3004  int32 vgroup_id = -1;
3005  int32 lone_vg_number = -1;
3006  int32 num_of_lones = 0;
3007  uint16 name_len = 0;
3008 
3009  file_id = Hopen (filename.c_str(), DFACC_READ, 0);
3010  if(file_id == FAIL)
3011  throw InternalErr(__FILE__,__LINE__,"Hopen failed.");
3012 
3013  status_n = Vstart (file_id);
3014  if(status_n == FAIL) {
3015  Hclose(file_id);
3016  throw InternalErr(__FILE__,__LINE__,"Vstart failed.");
3017  }
3018 
3019  string err_msg;
3020  bool unexpected_fail = false;
3021  //Get and print the names and class names of all the lone vgroups.
3022  // First, call Vlone with num_of_lones set to 0 to get the number of
3023  // lone vgroups in the file, but not to get their reference numbers.
3024  num_of_lones = Vlone (file_id, NULL, num_of_lones );
3025 
3026  //
3027  // Then, if there are any lone vgroups,
3028  if (num_of_lones > 0)
3029  {
3030  // Use the num_of_lones returned to allocate sufficient space for the
3031  // buffer ref_array to hold the reference numbers of all lone vgroups,
3032  vector<int32> ref_array;
3033  ref_array.resize(num_of_lones);
3034 
3035 
3036  // and call Vlone again to retrieve the reference numbers into
3037  // the buffer ref_array.
3038 
3039  num_of_lones = Vlone (file_id, &ref_array[0], num_of_lones);
3040 
3041  // Loop the name and class of each lone vgroup.
3042  for (lone_vg_number = 0; lone_vg_number < num_of_lones;
3043  lone_vg_number++)
3044  {
3045 
3046  // Attach to the current vgroup then get and display its
3047  // name and class. Note: the current vgroup must be detached before
3048  // moving to the next.
3049  vgroup_id = Vattach(file_id, ref_array[lone_vg_number], "r");
3050  if(vgroup_id == FAIL) {
3051  unexpected_fail = true;
3052  err_msg = string(ERR_LOC) + " Vattach failed. ";
3053  goto cleanFail;
3054  }
3055 
3056  status_32 = Vgetnamelen(vgroup_id, &name_len);
3057  if(status_32 == FAIL) {
3058  unexpected_fail = true;
3059  Vdetach(vgroup_id);
3060  err_msg = string(ERR_LOC) + " Vgetnamelen failed. ";
3061  goto cleanFail;
3062  }
3063 
3064  vector<char> vgroup_name;
3065  vgroup_name.resize(name_len+1);
3066  status_32 = Vgetname (vgroup_id, &vgroup_name[0]);
3067  if(status_32 == FAIL) {
3068  unexpected_fail = true;
3069  Vdetach(vgroup_id);
3070  err_msg = string(ERR_LOC) + " Vgetname failed. ";
3071  goto cleanFail;
3072  }
3073 
3074  status_32 = Vgetclassnamelen(vgroup_id, &name_len);
3075  if(status_32 == FAIL) {
3076  unexpected_fail = true;
3077  Vdetach(vgroup_id);
3078  err_msg = string(ERR_LOC) + " Vgetclassnamelen failed. ";
3079  goto cleanFail;
3080  }
3081 
3082  vector<char>vgroup_class;
3083  vgroup_class.resize(name_len+1);
3084  status_32 = Vgetclass (vgroup_id, &vgroup_class[0]);
3085  if(status_32 == FAIL) {
3086  unexpected_fail = true;
3087  Vdetach(vgroup_id);
3088  err_msg = string(ERR_LOC) + " Vgetclass failed. ";
3089  goto cleanFail;
3090  }
3091 
3092  string vgroup_name_str(vgroup_name.begin(),vgroup_name.end());
3093  vgroup_name_str = vgroup_name_str.substr(0,vgroup_name_str.size()-1);
3094 
3095  string vgroup_class_str(vgroup_class.begin(),vgroup_class.end());
3096  vgroup_class_str = vgroup_class_str.substr(0,vgroup_class_str.size()-1);
3097  try {
3098  if(vgroup_class_str =="GRID")
3099  map_eos2_one_object_attrs_wrapper(das,file_id,vgroup_id,vgroup_name_str,true);
3100  else if(vgroup_class_str =="SWATH")
3101  map_eos2_one_object_attrs_wrapper(das,file_id,vgroup_id,vgroup_name_str,false);
3102  }
3103  catch(...) {
3104  Vdetach(vgroup_id);
3105  Vend(file_id);
3106  Hclose(file_id);
3107  throw InternalErr(__FILE__,__LINE__,"map_eos2_one_object_attrs_wrapper failed.");
3108  }
3109  Vdetach (vgroup_id);
3110  }// for
3111  }// if
3112 
3113  //Terminate access to the V interface and close the file.
3114 cleanFail:
3115  Vend (file_id);
3116  Hclose (file_id);
3117  if(true == unexpected_fail)
3118  throw InternalErr(__FILE__,__LINE__,err_msg);
3119 
3120  return;
3121 
3122 }
3123 
3124 void HDFCFUtil::map_eos2_one_object_attrs_wrapper(libdap:: DAS &das,int32 file_id,int32 vgroup_id, const string& vgroup_name,bool is_grid) {
3125 
3126  //char attr_name[H4_MAX_NC_NAME]; //unused variable. SBL 2/7/20
3127 
3128  int32 num_gobjects = Vntagrefs (vgroup_id);
3129  if(num_gobjects < 0)
3130  throw InternalErr(__FILE__,__LINE__,"Cannot obtain the number of objects under a vgroup.");
3131 
3132  for(int i = 0; i<num_gobjects;i++) {
3133 
3134  int32 obj_tag, obj_ref;
3135  if (Vgettagref (vgroup_id, i, &obj_tag, &obj_ref) == FAIL)
3136  throw InternalErr(__FILE__,__LINE__,"Failed to obtain the tag and reference of an object under a vgroup.");
3137 
3138  if (Visvg (vgroup_id, obj_ref) == TRUE) {
3139 
3140  int32 object_attr_vgroup = Vattach(file_id,obj_ref,"r");
3141  if(object_attr_vgroup == FAIL)
3142  throw InternalErr(__FILE__,__LINE__,"Failed to attach an EOS2 vgroup.");
3143 
3144  uint16 name_len = 0;
3145  int32 status_32 = Vgetnamelen(object_attr_vgroup, &name_len);
3146  if(status_32 == FAIL) {
3147  Vdetach(object_attr_vgroup);
3148  throw InternalErr(__FILE__,__LINE__,"Failed to obtain an EOS2 vgroup name length.");
3149  }
3150  vector<char> attr_vgroup_name;
3151  attr_vgroup_name.resize(name_len+1);
3152  status_32 = Vgetname (object_attr_vgroup, &attr_vgroup_name[0]);
3153  if(status_32 == FAIL) {
3154  Vdetach(object_attr_vgroup);
3155  throw InternalErr(__FILE__,__LINE__,"Failed to obtain an EOS2 vgroup name. ");
3156  }
3157 
3158  string attr_vgroup_name_str(attr_vgroup_name.begin(),attr_vgroup_name.end());
3159  attr_vgroup_name_str = attr_vgroup_name_str.substr(0,attr_vgroup_name_str.size()-1);
3160 
3161  try {
3162  if(true == is_grid && attr_vgroup_name_str=="Grid Attributes"){
3163  map_eos2_one_object_attrs(das,file_id,object_attr_vgroup,vgroup_name);
3164  Vdetach(object_attr_vgroup);
3165  break;
3166  }
3167  else if(false == is_grid && attr_vgroup_name_str=="Swath Attributes") {
3168  map_eos2_one_object_attrs(das,file_id,object_attr_vgroup,vgroup_name);
3169  Vdetach(object_attr_vgroup);
3170  break;
3171  }
3172  }
3173  catch(...) {
3174  Vdetach(object_attr_vgroup);
3175  throw InternalErr(__FILE__,__LINE__,"Cannot map eos2 object attributes to DAP2.");
3176  }
3177  Vdetach(object_attr_vgroup);
3178 
3179  }
3180 
3181  }
3182 }
3183 
3184 void HDFCFUtil::map_eos2_one_object_attrs(libdap:: DAS &das,int32 file_id, int32 obj_attr_group_id, const string& vgroup_name) {
3185 
3186  int32 num_gobjects = Vntagrefs(obj_attr_group_id);
3187  if(num_gobjects < 0)
3188  throw InternalErr(__FILE__,__LINE__,"Cannot obtain the number of objects under a vgroup.");
3189 
3190  string vgroup_cf_name = HDFCFUtil::get_CF_string(vgroup_name);
3191  AttrTable *at = das.get_table(vgroup_cf_name);
3192  if(!at)
3193  at = das.add_table(vgroup_cf_name,new AttrTable);
3194 
3195  for(int i = 0; i<num_gobjects;i++) {
3196 
3197  int32 obj_tag, obj_ref;
3198  if (Vgettagref(obj_attr_group_id, i, &obj_tag, &obj_ref) == FAIL)
3199  throw InternalErr(__FILE__,__LINE__,"Failed to obtain the tag and reference of an object under a vgroup.");
3200 
3201  if(Visvs(obj_attr_group_id,obj_ref)) {
3202 
3203  int32 vdata_id = VSattach(file_id,obj_ref,"r");
3204  if(vdata_id == FAIL)
3205  throw InternalErr(__FILE__,__LINE__,"Failed to attach a vdata.");
3206 
3207  // EOS2 object vdatas are actually EOS2 object attributes.
3208  if(VSisattr(vdata_id)) {
3209 
3210  // According to EOS2 library, EOS2 number of field and record must be 1.
3211  int32 num_field = VFnfields(vdata_id);
3212  if(num_field !=1) {
3213  VSdetach(vdata_id);
3214  throw InternalErr(__FILE__,__LINE__,"Number of vdata field for an EOS2 object must be 1.");
3215  }
3216  int32 num_record = VSelts(vdata_id);
3217  if(num_record !=1){
3218  VSdetach(vdata_id);
3219  throw InternalErr(__FILE__,__LINE__,"Number of vdata record for an EOS2 object must be 1.");
3220  }
3221  char vdata_name[VSNAMELENMAX];
3222  if(VSQueryname(vdata_id,vdata_name) == FAIL) {
3223  VSdetach(vdata_id);
3224  throw InternalErr(__FILE__,__LINE__,"Failed to obtain EOS2 object vdata name.");
3225  }
3226  string vdatanamestr(vdata_name);
3227  string vdataname_cfstr = HDFCFUtil::get_CF_string(vdatanamestr);
3228  int32 fieldsize = VFfieldesize(vdata_id,0);
3229  if(fieldsize == FAIL) {
3230  VSdetach(vdata_id);
3231  throw InternalErr(__FILE__,__LINE__,"Failed to obtain EOS2 object vdata field size.");
3232  }
3233 
3234  char* fieldname = VFfieldname(vdata_id,0);
3235  if(fieldname == NULL) {
3236  VSdetach(vdata_id);
3237  throw InternalErr(__FILE__,__LINE__,"Failed to obtain EOS2 object vdata field name.");
3238  }
3239  int32 fieldtype = VFfieldtype(vdata_id,0);
3240  if(fieldtype == FAIL) {
3241  VSdetach(vdata_id);
3242  throw InternalErr(__FILE__,__LINE__,"Failed to obtain EOS2 object vdata field type.");
3243  }
3244 
3245  if(VSsetfields(vdata_id,fieldname) == FAIL) {
3246  VSdetach(vdata_id);
3247  throw InternalErr(__FILE__,__LINE__,"EOS2 object vdata: VSsetfields failed.");
3248  }
3249 
3250  vector<char> vdata_value;
3251  vdata_value.resize(fieldsize);
3252  if(VSread(vdata_id,(uint8*)&vdata_value[0],1,FULL_INTERLACE) == FAIL) {
3253  VSdetach(vdata_id);
3254  throw InternalErr(__FILE__,__LINE__,"EOS2 object vdata: VSread failed.");
3255  }
3256 
3257  // Map the attributes to DAP2.
3258  if(fieldtype == DFNT_UCHAR || fieldtype == DFNT_CHAR){
3259  string tempstring(vdata_value.begin(),vdata_value.end());
3260  // Remove the NULL term
3261  string tempstring2 = string(tempstring.c_str());
3262  at->append_attr(vdataname_cfstr,"String",HDFCFUtil::escattr(tempstring2));
3263  }
3264  else {
3265  string print_rep = HDFCFUtil::print_attr(fieldtype, 0, (void*) &vdata_value[0]);
3266  at->append_attr(vdataname_cfstr, HDFCFUtil::print_type(fieldtype), print_rep);
3267  }
3268 
3269  }
3270  VSdetach(vdata_id);
3271 
3272  }
3273  }
3274 
3275  return;
3276 }
3277 
3278 // Part of a large fix for attributes. Escaping the values of the attributes
3279 // may have been a bad idea. It breaks using JSON, for example. If this is a
3280 // bad idea - to turn of escaping - then we'll have to figure out how to store
3281 // 'serialized JSON' in attributes because it's being used in netcdf/hdf files.
3282 // If we stick with this, there's clearly a more performant solution - eliminate
3283 // the calls to this code.
3284 // jhrg 6/25/21
3285 #define ESCAPE_STRING_ATTRIBUTES 0
3286 
3287 string HDFCFUtil::escattr(string s)
3288 {
3289  const string printable = " ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789~`!@#$%^&*()_-+={[}]|\\:;<,>.?/'\"\n\t\r";
3290  const string ESC = "\\";
3291 #if ESCAPE_STRING_ATTRIBUTES
3292  const string DOUBLE_ESC = ESC + ESC;
3293  const string QUOTE = "\"";
3294  const string ESCQUOTE = ESC + QUOTE;
3295 
3296  // escape \ with a second backslash
3297  size_t ind = 0;
3298  while ((ind = s.find(ESC, ind)) != string::npos) {
3299  s.replace(ind, 1, DOUBLE_ESC);
3300  ind += DOUBLE_ESC.length();
3301  }
3302 
3303  // escape " with backslash
3304  ind = 0;
3305  while ((ind = s.find(QUOTE, ind)) != string::npos) {
3306  //comment out the following line since it wastes the CPU operation.
3307  s.replace(ind, 1, ESCQUOTE);
3308  ind += ESCQUOTE.length();
3309  }
3310 #endif
3311 
3312  size_t ind = 0;
3313  while ((ind = s.find_first_not_of(printable, ind)) != string::npos) {
3314  s.replace(ind, 1, ESC + octstring(s[ind]));
3315  }
3316 
3317  return s;
3318 }
3319 
3320 
3321 void HDFCFUtil::parser_trmm_v7_gridheader(const vector<char>& value,
3322  int& latsize, int&lonsize,
3323  float& lat_start, float& lon_start,
3324  float& lat_res, float& lon_res,
3325  bool check_reg_orig ){
3326 
3327  float lat_north = 0.;
3328  float lat_south = 0.;
3329  float lon_east = 0.;
3330  float lon_west = 0.;
3331 
3332  vector<string> ind_elems;
3333  char sep='\n';
3334  HDFCFUtil::Split(&value[0],sep,ind_elems);
3335  // The number of elements in the GridHeader is 9.
3336  //The string vector will add a leftover. So the size should be 10.
3337  // For the MacOS clang compiler, the string vector size may become 11.
3338  // So we change the condition to be "<9" is wrong.
3339  if(ind_elems.size()<9)
3340  throw InternalErr(__FILE__,__LINE__,"The number of elements in the TRMM level 3 GridHeader is not right.");
3341 
3342  if(false == check_reg_orig) {
3343  if (0 != ind_elems[1].find("Registration=CENTER"))
3344  throw InternalErr(__FILE__,__LINE__,"The TRMM grid registration is not center.");
3345  }
3346 
3347  if (0 == ind_elems[2].find("LatitudeResolution")){
3348 
3349  size_t equal_pos = ind_elems[2].find_first_of('=');
3350  if(string::npos == equal_pos)
3351  throw InternalErr(__FILE__,__LINE__,"Cannot find latitude resolution for TRMM level 3 products");
3352 
3353  size_t scolon_pos = ind_elems[2].find_first_of(';');
3354  if(string::npos == scolon_pos)
3355  throw InternalErr(__FILE__,__LINE__,"Cannot find latitude resolution for TRMM level 3 products");
3356  if (equal_pos < scolon_pos){
3357 
3358  string latres_str = ind_elems[2].substr(equal_pos+1,scolon_pos-equal_pos-1);
3359  lat_res = strtof(latres_str.c_str(),NULL);
3360  }
3361  else
3362  throw InternalErr(__FILE__,__LINE__,"latitude resolution is not right for TRMM level 3 products");
3363  }
3364  else
3365  throw InternalErr(__FILE__,__LINE__,"The TRMM grid LatitudeResolution doesn't exist.");
3366 
3367  if (0 == ind_elems[3].find("LongitudeResolution")){
3368 
3369  size_t equal_pos = ind_elems[3].find_first_of('=');
3370  if(string::npos == equal_pos)
3371  throw InternalErr(__FILE__,__LINE__,"Cannot find longitude resolution for TRMM level 3 products");
3372 
3373  size_t scolon_pos = ind_elems[3].find_first_of(';');
3374  if(string::npos == scolon_pos)
3375  throw InternalErr(__FILE__,__LINE__,"Cannot find longitude resolution for TRMM level 3 products");
3376  if (equal_pos < scolon_pos){
3377  string lonres_str = ind_elems[3].substr(equal_pos+1,scolon_pos-equal_pos-1);
3378  lon_res = strtof(lonres_str.c_str(),NULL);
3379  }
3380  else
3381  throw InternalErr(__FILE__,__LINE__,"longitude resolution is not right for TRMM level 3 products");
3382  }
3383  else
3384  throw InternalErr(__FILE__,__LINE__,"The TRMM grid LongitudeResolution doesn't exist.");
3385 
3386  if (0 == ind_elems[4].find("NorthBoundingCoordinate")){
3387 
3388  size_t equal_pos = ind_elems[4].find_first_of('=');
3389  if(string::npos == equal_pos)
3390  throw InternalErr(__FILE__,__LINE__,"Cannot find latitude resolution for TRMM level 3 products");
3391 
3392  size_t scolon_pos = ind_elems[4].find_first_of(';');
3393  if(string::npos == scolon_pos)
3394  throw InternalErr(__FILE__,__LINE__,"Cannot find latitude resolution for TRMM level 3 products");
3395  if (equal_pos < scolon_pos){
3396  string north_bounding_str = ind_elems[4].substr(equal_pos+1,scolon_pos-equal_pos-1);
3397  lat_north = strtof(north_bounding_str.c_str(),NULL);
3398  }
3399  else
3400  throw InternalErr(__FILE__,__LINE__,"NorthBoundingCoordinate is not right for TRMM level 3 products");
3401 
3402  }
3403  else
3404  throw InternalErr(__FILE__,__LINE__,"The TRMM grid NorthBoundingCoordinate doesn't exist.");
3405 
3406  if (0 == ind_elems[5].find("SouthBoundingCoordinate")){
3407 
3408  size_t equal_pos = ind_elems[5].find_first_of('=');
3409  if(string::npos == equal_pos)
3410  throw InternalErr(__FILE__,__LINE__,"Cannot find south bound coordinate for TRMM level 3 products");
3411 
3412  size_t scolon_pos = ind_elems[5].find_first_of(';');
3413  if(string::npos == scolon_pos)
3414  throw InternalErr(__FILE__,__LINE__,"Cannot find south bound coordinate for TRMM level 3 products");
3415  if (equal_pos < scolon_pos){
3416  string lat_south_str = ind_elems[5].substr(equal_pos+1,scolon_pos-equal_pos-1);
3417  lat_south = strtof(lat_south_str.c_str(),NULL);
3418  }
3419  else
3420  throw InternalErr(__FILE__,__LINE__,"south bound coordinate is not right for TRMM level 3 products");
3421 
3422  }
3423  else
3424  throw InternalErr(__FILE__,__LINE__,"The TRMM grid SouthBoundingCoordinate doesn't exist.");
3425 
3426  if (0 == ind_elems[6].find("EastBoundingCoordinate")){
3427 
3428  size_t equal_pos = ind_elems[6].find_first_of('=');
3429  if(string::npos == equal_pos)
3430  throw InternalErr(__FILE__,__LINE__,"Cannot find south bound coordinate for TRMM level 3 products");
3431 
3432  size_t scolon_pos = ind_elems[6].find_first_of(';');
3433  if(string::npos == scolon_pos)
3434  throw InternalErr(__FILE__,__LINE__,"Cannot find south bound coordinate for TRMM level 3 products");
3435  if (equal_pos < scolon_pos){
3436  string lon_east_str = ind_elems[6].substr(equal_pos+1,scolon_pos-equal_pos-1);
3437  lon_east = strtof(lon_east_str.c_str(),NULL);
3438  }
3439  else
3440  throw InternalErr(__FILE__,__LINE__,"south bound coordinate is not right for TRMM level 3 products");
3441 
3442  }
3443  else
3444  throw InternalErr(__FILE__,__LINE__,"The TRMM grid EastBoundingCoordinate doesn't exist.");
3445 
3446  if (0 == ind_elems[7].find("WestBoundingCoordinate")){
3447 
3448  size_t equal_pos = ind_elems[7].find_first_of('=');
3449  if(string::npos == equal_pos)
3450  throw InternalErr(__FILE__,__LINE__,"Cannot find south bound coordinate for TRMM level 3 products");
3451 
3452  size_t scolon_pos = ind_elems[7].find_first_of(';');
3453  if(string::npos == scolon_pos)
3454  throw InternalErr(__FILE__,__LINE__,"Cannot find south bound coordinate for TRMM level 3 products");
3455  if (equal_pos < scolon_pos){
3456  string lon_west_str = ind_elems[7].substr(equal_pos+1,scolon_pos-equal_pos-1);
3457  lon_west = strtof(lon_west_str.c_str(),NULL);
3458  }
3459  else
3460  throw InternalErr(__FILE__,__LINE__,"south bound coordinate is not right for TRMM level 3 products");
3461 
3462  }
3463  else
3464  throw InternalErr(__FILE__,__LINE__,"The TRMM grid WestBoundingCoordinate doesn't exist.");
3465 
3466  if (false == check_reg_orig) {
3467  if (0 != ind_elems[8].find("Origin=SOUTHWEST"))
3468  throw InternalErr(__FILE__,__LINE__,"The TRMM grid origin is not SOUTHWEST.");
3469  }
3470 
3471  // Since we only treat the case when the Registration is center, so the size should be the
3472  // regular number size - 1.
3473  latsize =(int)((lat_north-lat_south)/lat_res);
3474  lonsize =(int)((lon_east-lon_west)/lon_res);
3475  lat_start = lat_south;
3476  lon_start = lon_west;
3477 }
3478 
3479 // Somehow the conversion of double to c++ string with sprintf causes the memory error in
3480 // the testing code. I used the following code to do the conversion. Most part of the code
3481 // in reverse, int_to_str and dtoa are adapted from geeksforgeeks.org
3482 
3483 // reverses a string 'str' of length 'len'
3484 void HDFCFUtil::rev_str(char *str, int len)
3485 {
3486  int i=0;
3487  int j=len-1;
3488  int temp = 0;
3489  while (i<j)
3490  {
3491  temp = str[i];
3492  str[i] = str[j];
3493  str[j] = temp;
3494  i++;
3495  j--;
3496  }
3497 }
3498 
3499 // Converts a given integer x to string str[]. d is the number
3500 // of digits required in output. If d is more than the number
3501 // of digits in x, then 0s are added at the beginning.
3502 int HDFCFUtil::int_to_str(int x, char str[], int d)
3503 {
3504  int i = 0;
3505  while (x)
3506  {
3507  str[i++] = (x%10) + '0';
3508  x = x/10;
3509  }
3510 
3511  // If number of digits required is more, then
3512  // add 0s at the beginning
3513  while (i < d)
3514  str[i++] = '0';
3515 
3516  rev_str(str, i);
3517  str[i] = '\0';
3518  return i;
3519 }
3520 
3521 // Converts a double floating point number to string.
3522 void HDFCFUtil::dtoa(double n, char *res, int afterpoint)
3523 {
3524  // Extract integer part
3525  int ipart = (int)n;
3526 
3527  // Extract the double part
3528  double fpart = n - (double)ipart;
3529 
3530  // convert integer part to string
3531  int i = int_to_str(ipart, res, 0);
3532 
3533  // check for display option after point
3534  if (afterpoint != 0)
3535  {
3536  res[i] = '.'; // add dot
3537 
3538  // Get the value of fraction part upto given no.
3539  // of points after dot. The third parameter is needed
3540  // to handle cases like 233.007
3541  fpart = fpart * pow(10, afterpoint);
3542 
3543  // A round-error of 1 is found when casting to the integer for some numbers.
3544  // We need to correct it.
3545  int final_fpart = (int)fpart;
3546  if(fpart -(int)fpart >0.5)
3547  final_fpart = (int)fpart +1;
3548  int_to_str(final_fpart, res + i + 1, afterpoint);
3549  }
3550 }
3551 
3552 
3553 string HDFCFUtil::get_double_str(double x,int total_digit,int after_point) {
3554 
3555  string str;
3556  if(x!=0) {
3557  vector<char> res;
3558  res.resize(total_digit);
3559  for(int i = 0; i<total_digit;i++)
3560  res[i] = '\0';
3561  if (x<0) {
3562  str.push_back('-');
3563  dtoa(-x,&res[0],after_point);
3564  for(int i = 0; i<total_digit;i++) {
3565  if(res[i] != '\0')
3566  str.push_back(res[i]);
3567  }
3568  }
3569  else {
3570  dtoa(x, &res[0], after_point);
3571  for(int i = 0; i<total_digit;i++) {
3572  if(res[i] != '\0')
3573  str.push_back(res[i]);
3574  }
3575  }
3576 
3577  }
3578  else
3579  str.push_back('0');
3580 
3581  return str;
3582 
3583 }
3584 
3585 string HDFCFUtil::get_int_str(int x) {
3586 
3587  string str;
3588  if(x > 0 && x <10)
3589  str.push_back(x+'0');
3590 
3591  else if (x >10 && x<100) {
3592  str.push_back(x/10+'0');
3593  str.push_back(x%10+'0');
3594  }
3595  else {
3596  int num_digit = 0;
3597  int abs_x = (x<0)?-x:x;
3598  while(abs_x/=10)
3599  num_digit++;
3600  if(x<=0)
3601  num_digit++;
3602  vector<char> buf;
3603  buf.resize(num_digit);
3604  snprintf(&buf[0],num_digit,"%d",x);
3605  str.assign(&buf[0]);
3606 
3607  }
3608 
3609  return str;
3610 
3611 }
3612 
3613 #if 0
3614 template<typename T>
3615 size_t HDFCFUtil::write_vector_to_file(const string & fname, const vector<T> &val, size_t dtypesize) {
3616 #endif
3617 #if 0
3618 size_t HDFCFUtil::write_vector_to_file(const string & fname, const vector<double> &val, size_t dtypesize) {
3619 
3620  size_t ret_val;
3621  FILE* pFile;
3622 cerr<<"Open a file with the name "<<fname<<endl;
3623  pFile = fopen(fname.c_str(),"wb");
3624  ret_val = fwrite(&val[0],dtypesize,val.size(),pFile);
3625 cerr<<"ret_val for write is "<<ret_val <<endl;
3626 //for (int i=0;i<val.size();i++)
3627 //cerr<<"val["<<i<<"] is "<<val[i]<<endl;
3628  fclose(pFile);
3629  return ret_val;
3630 }
3631 ssize_t HDFCFUtil::write_vector_to_file2(const string & fname, const vector<double> &val, size_t dtypesize) {
3632 
3633  ssize_t ret_val;
3634  //int fd = open(fname.c_str(),O_RDWR|O_CREAT|O_TRUNC,0666);
3635  int fd = open(fname.c_str(),O_RDWR|O_CREAT|O_EXCL,0666);
3636 cerr<<"The first val is "<<val[0] <<endl;
3637  ret_val = write(fd,&val[0],dtypesize*val.size());
3638  close(fd);
3639 cerr<<"ret_val for write is "<<ret_val <<endl;
3640 //for (int i=0;i<val.size();i++)
3641 //cerr<<"val["<<i<<"] is "<<val[i]<<endl;
3642  return ret_val;
3643 }
3644 #endif
3645 ssize_t HDFCFUtil::read_vector_from_file(int fd, vector<double> &val, size_t dtypesize) {
3646 
3647  ssize_t ret_val;
3648  ret_val = read(fd,&val[0],val.size()*dtypesize);
3649 
3650 #if 0
3651 cerr<<"Open a file with the name "<<fname<<endl;
3652  pFile = fopen(fname.c_str(),"wb");
3653  ret_val = fwrite(&val[0],dtypesize,val.size(),pFile);
3654 cerr<<"ret_val for write is "<<ret_val <<endl;
3655 //for (int i=0;i<val.size();i++)
3656 //cerr<<"val["<<i<<"] is "<<val[i]<<endl;
3657  fclose(pFile);
3658 #endif
3659  return ret_val;
3660 }
3661 // Need to wrap a 'read buffer' from a pure file call here since read() is also a DAP function to read DAP data.
3662 ssize_t HDFCFUtil::read_buffer_from_file(int fd, void*buf, size_t total_read) {
3663 
3664  ssize_t ret_val;
3665  ret_val = read(fd,buf,total_read);
3666 
3667  return ret_val;
3668 }
3669 void HDFCFUtil::close_fileid(int32 sdfd, int32 fileid,int32 gridfd, int32 swathfd,bool pass_fileid) {
3670 
3671  if(false == pass_fileid) {
3672  if(sdfd != -1)
3673  SDend(sdfd);
3674  if(fileid != -1)
3675  Hclose(fileid);
3676 #ifdef USE_HDFEOS2_LIB
3677  if(gridfd != -1)
3678  GDclose(gridfd);
3679  if(swathfd != -1)
3680  SWclose(swathfd);
3681 
3682 #endif
3683  }
3684 
3685 }
3686 
3687 // Obtain the cache name. Since AIRS version 6 level 3 all share the same latitude and longitude,
3688 // we provide one set of latitude and longitude cache files for all AIRS level 3 version 6 products.
3689 string HDFCFUtil::obtain_cache_fname(const string & fprefix, const string &fname, const string &vname) {
3690 
3691  string cache_fname = fprefix;
3692  string base_file_name = basename(fname);
3693  // Identify this file from product name: AIRS, product level: .L3. or .L2. and version .v6.
3694  if((base_file_name.size() >12)
3695  && (base_file_name.compare(0,4,"AIRS") == 0)
3696  && (base_file_name.find(".L3.")!=string::npos)
3697  && (base_file_name.find(".v6.")!=string::npos)
3698  && ((vname == "Latitude") ||(vname == "Longitude")))
3699  cache_fname = cache_fname +"AIRS"+".L3."+".v6."+vname;
3700  else
3701  cache_fname = cache_fname + base_file_name +"_"+vname;
3702 
3703  return cache_fname;
3704 }
3705 
3706 // The current DDS cache is only for some products of which object information
3707 // can be retrieved via HDF4 SDS APIs. Currently only AIRS version 6 products are supported.
3708 size_t HDFCFUtil::obtain_dds_cache_size(HDFSP::File*spf) {
3709 
3710  size_t total_bytes_written = 0;
3711  const vector<HDFSP::SDField *>& spsds = spf->getSD()->getFields();
3712  vector<HDFSP::SDField *>::const_iterator it_g;
3713  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
3714 
3715  // We will not handle when the SDS datatype is DFNT_CHAR now.
3716  if(DFNT_CHAR == (*it_g)->getType()) {
3717  total_bytes_written = 0;
3718  break;
3719  }
3720  else {
3721  // We need to store dimension names and variable names.
3722  int temp_rank = (*it_g)->getRank();
3723  const vector<HDFSP::Dimension*>& dims= (*it_g)->getDimensions();
3724  vector<HDFSP::Dimension*>::const_iterator it_d;
3725  for(it_d = dims.begin(); it_d != dims.end(); ++it_d)
3726  total_bytes_written += ((*it_d)->getName()).size()+1;
3727 
3728  total_bytes_written +=((*it_g)->getName()).size()+1;
3729 
3730  // Many a time variable new name is the same as variable name,
3731  // so we may just save one byte('\0') for such as a case.
3732  if(((*it_g)->getName()) != ((*it_g)->getNewName()))
3733  total_bytes_written +=((*it_g)->getNewName()).size()+1;
3734  else
3735  total_bytes_written +=1;
3736 
3737  // We need to store 4 integers: rank, variable datatype, SDS reference number, fieldtype
3738  total_bytes_written +=(temp_rank+4)*sizeof(int);
3739  }
3740  }
3741 
3742  if(total_bytes_written != 0)
3743  total_bytes_written +=1;
3744 
3745  return total_bytes_written;
3746 
3747 }
3748 
3749 // Write the DDS of the special SDS-only HDF to a cache.
3750 void HDFCFUtil::write_sp_sds_dds_cache(HDFSP::File* spf,FILE*dds_file,size_t total_bytes_dds_cache,const string &dds_filename) {
3751 
3752  BESDEBUG("h4"," Coming to write SDS DDS to a cache" << endl);
3753  char delimiter = '\0';
3754  char cend = '\n';
3755  size_t total_written_bytes_count = 0;
3756 
3757  // The buffer to hold information to write to a DDS cache file.
3758  vector<char>temp_buf;
3759  temp_buf.resize(total_bytes_dds_cache);
3760  char* temp_pointer = &temp_buf[0];
3761 
3762  const vector<HDFSP::SDField *>& spsds = spf->getSD()->getFields();
3763 
3764  //Read SDS
3765  vector<HDFSP::SDField *>::const_iterator it_g;
3766  for(it_g = spsds.begin(); it_g != spsds.end(); it_g++){
3767 
3768  // First, rank, fieldtype, SDS reference number, variable datatype, dimsize(rank)
3769  int sds_rank = (*it_g)->getRank();
3770  int sds_ref = (*it_g)->getFieldRef();
3771  int sds_dtype = (*it_g)->getType();
3772  int sds_ftype = (*it_g)->getFieldType();
3773 
3774  vector<int32>dimsizes;
3775  dimsizes.resize(sds_rank);
3776 
3777  // Size for each dimension
3778  const vector<HDFSP::Dimension*>& dims= (*it_g)->getDimensions();
3779  for(int i = 0; i < sds_rank; i++)
3780  dimsizes[i] = dims[i]->getSize();
3781 
3782  memcpy((void*)temp_pointer,(void*)&sds_rank,sizeof(int));
3783  temp_pointer +=sizeof(int);
3784  memcpy((void*)temp_pointer,(void*)&sds_ref,sizeof(int));
3785  temp_pointer +=sizeof(int);
3786  memcpy((void*)temp_pointer,(void*)&sds_dtype,sizeof(int));
3787  temp_pointer +=sizeof(int);
3788  memcpy((void*)temp_pointer,(void*)&sds_ftype,sizeof(int));
3789  temp_pointer +=sizeof(int);
3790 
3791  memcpy((void*)temp_pointer,(void*)&dimsizes[0],sds_rank*sizeof(int));
3792  temp_pointer +=sds_rank*sizeof(int);
3793 
3794  // total written bytes so far
3795  total_written_bytes_count +=(sds_rank+4)*sizeof(int);
3796 
3797  // Second, variable name,variable new name and SDS dim name(rank)
3798  // we need to a delimiter to distinguish each name.
3799  string temp_varname = (*it_g)->getName();
3800  vector<char>temp_val1(temp_varname.begin(),temp_varname.end());
3801  memcpy((void*)temp_pointer,(void*)&temp_val1[0],temp_varname.size());
3802  temp_pointer +=temp_varname.size();
3803  memcpy((void*)temp_pointer,&delimiter,1);
3804  temp_pointer +=1;
3805 
3806  total_written_bytes_count =total_written_bytes_count +(1+temp_varname.size());
3807 
3808  // To save the dds cache size and the speed to retrieve variable new name
3809  // we only save variable cf name when the variable cf name is not the
3810  // same as the variable name. When they are the same, a delimiter is saved for
3811  // variable cf name.
3812  if((*it_g)->getName() == (*it_g)->getNewName()){
3813  memcpy((void*)temp_pointer,&delimiter,1);
3814  temp_pointer +=1;
3815  total_written_bytes_count +=1;
3816  }
3817  else {
3818  string temp_varnewname = (*it_g)->getNewName();
3819  vector<char>temp_val2(temp_varnewname.begin(),temp_varnewname.end());
3820  memcpy((void*)temp_pointer,(void*)&temp_val2[0],temp_varnewname.size());
3821  temp_pointer +=temp_varnewname.size();
3822  memcpy((void*)temp_pointer,&delimiter,1);
3823  temp_pointer +=1;
3824  total_written_bytes_count =total_written_bytes_count +(1+temp_varnewname.size());
3825  }
3826 
3827  // Now the name for each dimensions.
3828  for(int i = 0; i < sds_rank; i++) {
3829  string temp_dimname = dims[i]->getName();
3830  vector<char>temp_val3(temp_dimname.begin(),temp_dimname.end());
3831  memcpy((void*)temp_pointer,(void*)&temp_val3[0],temp_dimname.size());
3832  temp_pointer +=temp_dimname.size();
3833  memcpy((void*)temp_pointer,&delimiter,1);
3834  temp_pointer +=1;
3835  total_written_bytes_count =total_written_bytes_count +(1+temp_dimname.size());
3836  }
3837  }
3838 
3839  memcpy((void*)temp_pointer,&cend,1);
3840  total_written_bytes_count +=1;
3841 
3842  if(total_written_bytes_count != total_bytes_dds_cache) {
3843  stringstream s_total_written_count;
3844  s_total_written_count << total_written_bytes_count;
3845  stringstream s_total_bytes_dds_cache;
3846  s_total_bytes_dds_cache << total_bytes_dds_cache;
3847  string msg = "DDs cached file "+ dds_filename +" buffer size should be " + s_total_bytes_dds_cache.str() ;
3848  msg = msg + ". But the real size written in the buffer is " + s_total_written_count.str();
3849  throw InternalErr (__FILE__, __LINE__,msg);
3850  }
3851 
3852  size_t bytes_really_written = fwrite((const void*)&temp_buf[0],1,total_bytes_dds_cache,dds_file);
3853 
3854  if(bytes_really_written != total_bytes_dds_cache) {
3855  stringstream s_expected_bytes;
3856  s_expected_bytes << total_bytes_dds_cache;
3857  stringstream s_really_written_bytes;
3858  s_really_written_bytes << bytes_really_written;
3859  string msg = "DDs cached file "+ dds_filename +" size should be " + s_expected_bytes.str() ;
3860  msg += ". But the real size written to the file is " + s_really_written_bytes.str();
3861  throw InternalErr (__FILE__, __LINE__,msg);
3862  }
3863 
3864 }
3865 
3866 // Read DDS of a special SDS-only HDF file from a cache.
3867 void HDFCFUtil::read_sp_sds_dds_cache(FILE* dds_file,libdap::DDS * dds_ptr,
3868  const std::string &cache_filename, const std::string &hdf4_filename) {
3869 
3870  BESDEBUG("h4"," Coming to read SDS DDS from a cache" << endl);
3871 
3872  // Check the cache file size.
3873  struct stat sb;
3874  if(stat(cache_filename.c_str(),&sb)!=0) {
3875  string err_mesg="The DDS cache file " + cache_filename;
3876  err_mesg = err_mesg + " doesn't exist. ";
3877  throw InternalErr(__FILE__,__LINE__,err_mesg);
3878  }
3879 
3880  size_t bytes_expected_read = (size_t)sb.st_size;
3881 
3882  // Allocate the buffer size based on the file size.
3883  vector<char> temp_buf;
3884  temp_buf.resize(bytes_expected_read);
3885  size_t bytes_really_read = fread((void*)&temp_buf[0],1,bytes_expected_read,dds_file);
3886 
3887  // Now bytes_really_read should be the same as bytes_expected_read if the element size is 1.
3888  if(bytes_really_read != bytes_expected_read) {
3889  stringstream s_bytes_really_read;
3890  s_bytes_really_read << bytes_really_read ;
3891  stringstream s_bytes_expected_read;
3892  s_bytes_expected_read << bytes_expected_read;
3893  string msg = "The expected bytes to read from DDS cache file " + cache_filename +" is " + s_bytes_expected_read.str();
3894  msg = msg + ". But the real read size from the buffer is " + s_bytes_really_read.str();
3895  throw InternalErr (__FILE__, __LINE__,msg);
3896  }
3897  char* temp_pointer = &temp_buf[0];
3898 
3899  char delimiter = '\0';
3900  // The end of the whole string.
3901  char cend = '\n';
3902  bool end_file_flag = false;
3903 
3904 
3905  do {
3906  int sds_rank = *((int *)(temp_pointer));
3907  temp_pointer = temp_pointer + sizeof(int);
3908 
3909  int sds_ref = *((int *)(temp_pointer));
3910  temp_pointer = temp_pointer + sizeof(int);
3911 
3912  int sds_dtype = *((int *)(temp_pointer));
3913  temp_pointer = temp_pointer + sizeof(int);
3914 
3915  int sds_ftype = *((int *)(temp_pointer));
3916  temp_pointer = temp_pointer + sizeof(int);
3917 
3918  vector<int32>dimsizes;
3919  if(sds_rank <=0)
3920  throw InternalErr (__FILE__, __LINE__,"SDS rank must be >0");
3921 
3922  dimsizes.resize(sds_rank);
3923  for (int i = 0; i <sds_rank;i++) {
3924  dimsizes[i] = *((int *)(temp_pointer));
3925  temp_pointer = temp_pointer + sizeof(int);
3926  }
3927 
3928  vector<string>dimnames;
3929  dimnames.resize(sds_rank);
3930  string varname,varnewname;
3931  for (int i = 0; i <sds_rank+2;i++) {
3932  vector<char> temp_vchar;
3933  char temp_char = *temp_pointer;
3934 
3935  // Only apply when varnewname is stored as the delimiter.
3936  if(temp_char == delimiter)
3937  temp_vchar.push_back(temp_char);
3938  while(temp_char !=delimiter) {
3939  temp_vchar.push_back(temp_char);
3940  temp_pointer++;
3941  temp_char = *temp_pointer;
3942  //temp_char = *(++temp_pointer); work
3943  //temp_char = *(temp_pointer++); not working
3944  }
3945  string temp_string(temp_vchar.begin(),temp_vchar.end());
3946  if(i == 0)
3947  varname = temp_string;
3948  else if( i == 1)
3949  varnewname = temp_string;
3950  else
3951  dimnames[i-2] = temp_string;
3952  temp_pointer++;
3953  }
3954 
3955  // If varnewname is only the delimiter, varname and varnewname is the same.
3956  if(varnewname[0] == delimiter)
3957  varnewname = varname;
3958  // Assemble DDS.
3959  // 1. Create a basetype
3960  BaseType *bt = NULL;
3961  switch(sds_dtype) {
3962 #define HANDLE_CASE(tid, type) \
3963  case tid: \
3964  bt = new (type)(varnewname,hdf4_filename); \
3965  break;
3966  HANDLE_CASE(DFNT_FLOAT32, HDFFloat32);
3967  HANDLE_CASE(DFNT_FLOAT64, HDFFloat64);
3968  HANDLE_CASE(DFNT_CHAR, HDFStr);
3969 #ifndef SIGNED_BYTE_TO_INT32
3970  HANDLE_CASE(DFNT_INT8, HDFByte);
3971 #else
3972  HANDLE_CASE(DFNT_INT8,HDFInt32);
3973 #endif
3974  HANDLE_CASE(DFNT_UINT8, HDFByte);
3975  HANDLE_CASE(DFNT_INT16, HDFInt16);
3976  HANDLE_CASE(DFNT_UINT16, HDFUInt16);
3977  HANDLE_CASE(DFNT_INT32, HDFInt32);
3978  HANDLE_CASE(DFNT_UINT32, HDFUInt32);
3979  HANDLE_CASE(DFNT_UCHAR8, HDFByte);
3980  default:
3981  throw InternalErr(__FILE__,__LINE__,"unsupported data type.");
3982 #undef HANDLE_CASE
3983  }
3984 
3985  if(NULL == bt)
3986  throw InternalErr(__FILE__,__LINE__,"Cannot create the basetype when creating DDS from a cache file.");
3987 
3988  SPType sptype = OTHERHDF;
3989 
3990  // sds_ftype indicates if this is a general data field or geolocation field.
3991  // 4 indicates this is a missing non-latlon geo-location fields.
3992  if(sds_ftype != 4){
3993  HDFSPArray_RealField *ar = NULL;
3994  try {
3995  // pass sds id as 0 since the sds id will be retrieved from SDStart if necessary.
3996  ar = new HDFSPArray_RealField(
3997  sds_rank,
3998  hdf4_filename,
3999  0,
4000  sds_ref,
4001  sds_dtype,
4002  sptype,
4003  varname,
4004  dimsizes,
4005  varnewname,
4006  bt);
4007  }
4008  catch(...) {
4009  delete bt;
4010  throw InternalErr(__FILE__,__LINE__,"Unable to allocate the HDFSPArray_RealField instance.");
4011  }
4012 
4013  for(int i = 0; i <sds_rank; i++)
4014  ar->append_dim(dimsizes[i],dimnames[i]);
4015  dds_ptr->add_var(ar);
4016  delete bt;
4017  delete ar;
4018  }
4019  else {
4020  HDFSPArrayMissGeoField *ar = NULL;
4021  if(sds_rank == 1) {
4022  try {
4023  ar = new HDFSPArrayMissGeoField(
4024  sds_rank,
4025  dimsizes[0],
4026  varnewname,
4027  bt);
4028  }
4029  catch(...) {
4030  delete bt;
4031  throw InternalErr(__FILE__,__LINE__,"Unable to allocate the HDFSPArray_RealField instance.");
4032  }
4033 
4034  ar->append_dim(dimsizes[0],dimnames[0]);
4035  dds_ptr->add_var(ar);
4036  delete bt;
4037  delete ar;
4038  }
4039  else
4040  throw InternalErr(__FILE__,__LINE__,"SDS rank must be 1 for the missing coordinate.");
4041  }
4042 
4043  if(*temp_pointer == cend)
4044  end_file_flag = true;
4045  if((temp_pointer - &temp_buf[0]) > (int)bytes_expected_read) {
4046  string msg = cache_filename + " doesn't have the end-line character at the end. The file may be corrupted.";
4047  throw InternalErr (__FILE__, __LINE__,msg);
4048  }
4049  } while(false == end_file_flag);
4050 
4051  dds_ptr->set_dataset_name(basename(hdf4_filename));
4052 }
4053 
4054 
4055 #if 0
4056 void HDFCFUtil::close_fileid(int32 sdfd, int32 fileid,int32 gridfd, int32 swathfd) {
4057 
4058  if(sdfd != -1)
4059  SDend(sdfd);
4060  if(fileid != -1)
4061  Hclose(fileid);
4062  if(gridfd != -1)
4063  GDclose(gridfd);
4064  if(swathfd != -1)
4065  SWclose(swathfd);
4066 
4067 }
4068 
4069 void HDFCFUtil::reset_fileid(int& sdfd, int& fileid,int& gridfd, int& swathfd) {
4070 
4071  sdfd = -1;
4072  fileid = -1;
4073  gridfd = -1;
4074  swathfd = -1;
4075 
4076 }
4077 #endif
static std::string lowercase(const std::string &s)
Definition: BESUtil.cc:206
const std::vector< Attribute * > & getAttributes() const
Get the attributes of this field.
Definition: HDFSP.h:315
int32 getType() const
Get the data type of this field.
Definition: HDFSP.h:309
const std::string & getNewName() const
Get the CF name(special characters replaced by underscores) of this field.
Definition: HDFSP.h:297
bool Has_Dim_NoScale_Field() const
This file has a field that is a SDS dimension but no dimension scale.
Definition: HDFSP.h:756
SD * getSD() const
Public interface to Obtain SD.
Definition: HDFSP.h:771
const std::vector< VDATA * > & getVDATAs() const
Public interface to Obtain Vdata.
Definition: HDFSP.h:777
SPType getSPType() const
Obtain special HDF4 product type.
Definition: HDFSP.h:749
const std::string & getPath() const
Obtain the path of the file.
Definition: HDFSP.h:765
One instance of this class represents one SDS object.
Definition: HDFSP.h:346
This class retrieves all SDS objects and SD file attributes.
Definition: HDFSP.h:558
const std::vector< SDField * > & getFields() const
Redundant member function.
Definition: HDFSP.h:577
const std::vector< Attribute * > & getAttributes() const
Public interface to obtain the SD(file) attributes.
Definition: HDFSP.h:583
Definition: HDFStr.h:51
void get_value(const std::string &s, std::string &val, bool &found)
Retrieve the value of a given key, if set.
Definition: TheBESKeys.cc:340
static TheBESKeys * TheKeys()
Definition: TheBESKeys.cc:71
static short obtain_type_size(int32)
Obtain datatype size.
Definition: HDFCFUtil.cc:450
static bool insert_map(std::map< std::string, std::string > &m, std::string key, std::string val)
Definition: HDFCFUtil.cc:148
static std::string print_attr(int32, int, void *)
Print attribute values in string.
Definition: HDFCFUtil.cc:268
static std::string print_type(int32)
Print datatype in string.
Definition: HDFCFUtil.cc:389
static void gen_unique_name(std::string &str, std::set< std::string > &namelist, int &clash_index)
Obtain the unique name for the clashed names and save it to set namelist.
Definition: HDFCFUtil.cc:194
static void Handle_NameClashing(std::vector< std::string > &newobjnamelist)
General routines to handle name clashings.
Definition: HDFCFUtil.cc:260
static void correct_scale_offset_type(libdap::AttrTable *at)
Definition: HDFCFUtil.cc:614
static std::string get_CF_string(std::string s)
Change special characters to "_".
Definition: HDFCFUtil.cc:164
static void Split(const char *s, int len, char sep, std::vector< std::string > &names)
Definition: HDFCFUtil.cc:80
static std::string escattr(std::string s)
Definition: HDFCFUtil.cc:3287
static void correct_fvalue_type(libdap::AttrTable *at, int32 dtype)
Definition: HDFCFUtil.cc:547
static void close_fileid(int32 sdfd, int32 file_id, int32 gridfd, int32 swathfd, bool pass_fileid_key)
Definition: HDFCFUtil.cc:3669
Definition: HDFCFUtil.h:52