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