bes  Updated for version 3.20.10
HDF5Array.cc
Go to the documentation of this file.
1 // This file is part of hdf5_handler a HDF5 file handler for the OPeNDAP
2 // data server.
3 
4 // Author: Hyo-Kyung Lee <hyoklee@hdfgroup.org> and Muqun Yang
5 // <myang6@hdfgroup.org>
6 
7 // Copyright (c) 2009-2016 The HDF Group, Inc. and OPeNDAP, Inc.
8 //
9 // This is H5free_memory software; you can redistribute it and/or modify it under the
10 // terms of the GNU Lesser General Public License as published by the Free
11 // Software Foundation; either version 2.1 of the License, or (at your
12 // option) any later version.
13 //
14 // This software is distributed in the hope that it will be useful, but
15 // WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
16 // or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17 // License for more details.
18 //
19 // You should have received a copy of the GNU Lesser General Public
20 // License along with this library; if not, write to the Free Software
21 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
22 //
23 // You can contact OPeNDAP, Inc. at PO Box 112, Saunderstown, RI. 02874-0112.
24 // You can contact The HDF Group, Inc. at 1800 South Oak Street,
25 // Suite 203, Champaign, IL 61820
33 
34 
35 #include "config_hdf5.h"
36 
37 #include <iostream>
38 #include <memory>
39 #include <sstream>
40 #include <algorithm>
41 #include <ctype.h>
42 
43 #include <BESDebug.h>
44 #include <libdap/Error.h>
45 #include <libdap/InternalErr.h>
46 
47 #include "HDF5Array.h"
48 #include "HDF5Structure.h"
49 #include "HDF5Str.h"
50 
51 using namespace std;
52 using namespace libdap;
53 
55  return new HDF5Array(*this);
56 }
57 
58 HDF5Array::HDF5Array(const string & n, const string &d, BaseType * v) :
59  Array(n, d, v) {
60  d_num_dim = 0;
61  d_num_elm = 0;
62  d_memneed = 0;
63 }
64 
65 HDF5Array::~HDF5Array() {
66 }
67 
68 int HDF5Array::format_constraint(int *offset, int *step, int *count) {
69 
70  // For the 0-length array case, just return 0.
71  if(length() == 0)
72  return 0;
73 
74  //long nels = 1;
75  int nels = 1;
76  int id = 0;
77 
78  Dim_iter p = dim_begin();
79 
80  while (p != dim_end()) {
81 
82  int start = dimension_start(p, true);
83  int stride = dimension_stride(p, true);
84  int stop = dimension_stop(p, true);
85 
86  // Check for empty constraint
87  if (start > stop) {
88  ostringstream oss;
89 
90  oss << "Array/Grid hyperslab start point "<< start <<
91  " is greater than stop point " << stop <<".";
92  throw Error(malformed_expr, oss.str());
93  }
94 
95  offset[id] = start;
96  step[id] = stride;
97  count[id] = ((stop - start) / stride) + 1; // count of elements
98  nels *= count[id]; // total number of values for variable
99 
100  BESDEBUG("h5",
101  "=format_constraint():"
102  << "id=" << id << " offset=" << offset[id]
103  << " step=" << step[id]
104  << " count=" << count[id]
105  << endl);
106 
107  id++;
108  p++;
109  }
110 
111  return nels;
112 }
113 
114 
116 {
117  BESDEBUG("h5",
118  ">read() dataset=" << dataset()
119  << " dimension=" << d_num_dim
120  << " data_size=" << d_memneed << " length=" << length()
121  << endl);
122 
123  hid_t file_id = H5Fopen(dataset().c_str(),H5F_ACC_RDONLY,H5P_DEFAULT);
124 
125  BESDEBUG("h5","after H5Fopen "<<endl);
126  BESDEBUG("h5","variable name is "<<name() <<endl);
127  BESDEBUG("h5","variable path is "<<var_path <<endl);
128 
129  hid_t dset_id = -1;
130 
131  if(true == is_dap4())
132  dset_id = H5Dopen2(file_id,var_path.c_str(),H5P_DEFAULT);
133  else
134  dset_id = H5Dopen2(file_id,name().c_str(),H5P_DEFAULT);
135 
136  BESDEBUG("h5","after H5Dopen2 "<<endl);
137  // Leave the following code. We may replace the struct DS and DSattr(see hdf5_handler.h)
138 #if 0
139  hid_t dspace_id = H5Dget_space(dset_id);
140  if(dspace_id < 0) {
141  H5Dclose(dset_id);
142  H5Fclose(file_id);
143  throw InternalErr(__FILE__,__LINE__, "Fail to obtain the dataspace .");
144  }
145 
146  int num_dim = H5Sget_simple_extent_ndims(dspace_id);
147  if(num_dim < 0) {
148  H5Sclose(dspace_id);
149  H5Dclose(dset_id);
150  H5Fclose(file_id);
151  throw InternalErr(__FILE__,__LINE__, "Fail to obtain the datatype .");
152  }
153 
154  H5Sclose(dspace_id);
155 #endif
156 
157  hid_t dtype_id = H5Dget_type(dset_id);
158  if(dtype_id < 0) {
159  H5Dclose(dset_id);
160  H5Fclose(file_id);
161  throw InternalErr(__FILE__,__LINE__, "Fail to obtain the datatype .");
162  }
163 
164 
165  vector<int> offset(d_num_dim);
166  vector<int> count(d_num_dim);
167  vector<int> step(d_num_dim);
168  int nelms = format_constraint(&offset[0], &step[0], &count[0]); // Throws Error.
169  vector<char>values;
170 
171 
172  // We only map the reference to URL when the dataset is an array of reference.
173  if (get_dap_type(dtype_id,is_dap4()) == "Url") {
174  bool ret_ref = false;
175  try {
176  ret_ref = m_array_of_reference(dset_id,dtype_id);
177  H5Tclose(dtype_id);
178  H5Dclose(dset_id);
179  H5Fclose(file_id);
180 
181  }
182  catch(...) {
183  H5Tclose(dtype_id);
184  H5Dclose(dset_id);
185  H5Fclose(file_id);
186  throw;
187 
188  }
189  return ret_ref;
190  }
191 
192  try {
193  do_array_read(dset_id,dtype_id,values,false,0,nelms,&offset[0],&count[0],&step[0]);
194  }
195  catch(...) {
196  H5Tclose(dtype_id);
197  H5Dclose(dset_id);
198  H5Fclose(file_id);
199  throw;
200  }
201 
202  H5Tclose(dtype_id);
203  H5Dclose(dset_id);
204  H5Fclose(file_id);
205 
206  return true;
207 }
208 
209 void HDF5Array::do_array_read(hid_t dset_id,hid_t dtype_id,vector<char>&values,bool has_values,int values_offset,
210  int nelms,int* offset,int* count, int* step)
211 {
212 
213  H5T_class_t tcls = H5Tget_class(dtype_id);
214 
215  if(H5T_COMPOUND == tcls)
216  m_array_of_structure(dset_id,values,has_values,values_offset,nelms,offset,count,step);
217  else if(H5T_INTEGER == tcls || H5T_FLOAT == tcls || H5T_STRING == tcls)
218  m_array_of_atomic(dset_id,dtype_id,nelms,offset,count,step);
219  else {
220  throw InternalErr(__FILE__,__LINE__,"Fail to read the data for Unsupported datatype.");
221  }
222 
223 }
224 
225 void HDF5Array:: m_array_of_atomic(hid_t dset_id, hid_t dtype_id,
226  int nelms,int* offset,int* count, int* step)
227 {
228 
229  hid_t memtype = -1;
230  if((memtype = H5Tget_native_type(dtype_id, H5T_DIR_ASCEND))<0) {
231  throw InternalErr (__FILE__, __LINE__, "Fail to obtain memory datatype.");
232  }
233 
234  // First handle variable-length string
235  if (H5Tis_variable_str(memtype) && H5Tget_class(memtype) == H5T_STRING) {
236 
237  vector<hsize_t> hoffset;
238  vector<hsize_t>hcount;
239  vector<hsize_t>hstep;
240  hoffset.resize(d_num_dim);
241  hcount.resize(d_num_dim);
242  hstep.resize(d_num_dim);
243  for (int i = 0; i <d_num_dim; i++) {
244  hoffset[i] = (hsize_t) offset[i];
245  hcount[i] = (hsize_t) count[i];
246  hstep[i] = (hsize_t) step[i];
247  }
248 
249  vector<string>finstrval;
250  finstrval.resize(nelms);
251  try {
252  read_vlen_string(dset_id, nelms, &hoffset[0], &hstep[0], &hcount[0],finstrval);
253  }
254  catch(...) {
255  H5Tclose(memtype);
256  throw InternalErr(__FILE__,__LINE__,"Fail to read variable-length string.");
257  }
258  set_value(finstrval,nelms);
259  H5Tclose(memtype);
260  return ;
261  }
262 
263  try {
264  if (nelms == d_num_elm) {
265 
266  vector<char> convbuf(d_memneed);
267  get_data(dset_id, (void *) &convbuf[0]);
268 
269  // Check if a Signed Byte to Int16 conversion is necessary, this is only valid for DAP2.
270  if(false == is_dap4()) {
271  if (1 == H5Tget_size(memtype) && H5T_SGN_2 == H5Tget_sign(memtype))
272  {
273  vector<short> convbuf2(nelms);
274  for (int i = 0; i < nelms; i++) {
275  convbuf2[i] = (signed char) (convbuf[i]);
276  BESDEBUG("h5", "convbuf[" << i << "]="
277  << (signed char)convbuf[i] << endl);
278  BESDEBUG("h5", "convbuf2[" << i << "]="
279  << convbuf2[i] << endl)
280  ;
281  }
282  // Libdap will generate the wrong output.
283  m_intern_plain_array_data((char*) &convbuf2[0],memtype);
284  }
285  else
286  m_intern_plain_array_data(&convbuf[0],memtype);
287  }
288  else
289  m_intern_plain_array_data(&convbuf[0],memtype);
290  } // end of "if (nelms == d_num_elm)"
291  else {
292  size_t data_size = nelms * H5Tget_size(memtype);
293  if (data_size == 0) {
294  throw InternalErr(__FILE__, __LINE__, "get_size failed");
295  }
296  vector<char> convbuf(data_size);
297  get_slabdata(dset_id, &offset[0], &step[0], &count[0], d_num_dim, &convbuf[0]);
298 
299  // Check if a Signed Byte to Int16 conversion is necessary.
300  if(false == is_dap4()){
301  if (1 == H5Tget_size(memtype) && H5T_SGN_2 == H5Tget_sign(memtype)) {
302  vector<short> convbuf2(data_size);
303  for (int i = 0; i < (int)data_size; i++) {
304  convbuf2[i] = static_cast<signed char> (convbuf[i]);
305  }
306  m_intern_plain_array_data((char*) &convbuf2[0],memtype);
307  }
308  else {
309  m_intern_plain_array_data(&convbuf[0],memtype);
310  }
311  }
312  else
313  m_intern_plain_array_data(&convbuf[0],memtype);
314 
315  }
316  H5Tclose(memtype);
317  }
318  catch (...) {
319  H5Tclose(memtype);
320  throw;
321  }
322 
323 }
324 
325 bool HDF5Array::m_array_of_structure(hid_t dsetid, vector<char>&values,bool has_values,int values_offset,
326  int nelms,int* offset,int* count, int* step) {
327 
328  BESDEBUG("h5", "=read() Array of Structure length=" << length() << endl);
329 
330  hid_t mspace = -1;
331  hid_t memtype = -1;
332  hid_t dtypeid = -1;
333  size_t ty_size = -1;
334 
335  if((dtypeid = H5Dget_type(dsetid)) < 0)
336  throw InternalErr (__FILE__, __LINE__, "Cannot obtain the datatype.");
337 
338  if((memtype = H5Tget_native_type(dtypeid, H5T_DIR_ASCEND))<0) {
339  H5Tclose(dtypeid);
340  throw InternalErr (__FILE__, __LINE__, "Fail to obtain memory datatype.");
341  }
342 
343  ty_size = H5Tget_size(memtype);
344  if (ty_size == 0) {
345  H5Tclose(memtype);
346  H5Tclose(dtypeid);
347  throw InternalErr (__FILE__, __LINE__,"Fail to obtain the size of HDF5 compound datatype.");
348  }
349 
350  if(false == has_values) {
351 
352  hid_t dspace = -1;
353 
354  if ((dspace = H5Dget_space(dsetid))<0) {
355  H5Tclose(memtype);
356  H5Tclose(dtypeid);
357  throw InternalErr (__FILE__, __LINE__, "Cannot obtain data space.");
358  }
359 
360  d_num_dim = H5Sget_simple_extent_ndims(dspace);
361  if(d_num_dim < 0) {
362  H5Tclose(memtype);
363  H5Tclose(dtypeid);
364  H5Sclose(dspace);
365  throw InternalErr (__FILE__, __LINE__, "Cannot obtain the number of dimensions of the data space.");
366  }
367 
368  vector<hsize_t> hoffset;
369  vector<hsize_t>hcount;
370  vector<hsize_t>hstep;
371  hoffset.resize(d_num_dim);
372  hcount.resize(d_num_dim);
373  hstep.resize(d_num_dim);
374  for (int i = 0; i <d_num_dim; i++) {
375  hoffset[i] = (hsize_t) offset[i];
376  hcount[i] = (hsize_t) count[i];
377  hstep[i] = (hsize_t) step[i];
378  }
379 
380  if (H5Sselect_hyperslab(dspace, H5S_SELECT_SET,
381  &hoffset[0], &hstep[0],
382  &hcount[0], NULL) < 0) {
383  H5Tclose(memtype);
384  H5Tclose(dtypeid);
385  H5Sclose(dspace);
386  throw InternalErr (__FILE__, __LINE__, "Cannot generate the hyperslab of the HDF5 dataset.");
387  }
388 
389  mspace = H5Screate_simple(d_num_dim, &hcount[0],NULL);
390  if (mspace < 0) {
391  H5Sclose(dspace);
392  H5Tclose(memtype);
393  H5Tclose(dtypeid);
394  throw InternalErr (__FILE__, __LINE__, "Cannot create the memory space.");
395  }
396 
397  values.resize(nelms*ty_size);
398  hid_t read_ret = -1;
399  read_ret = H5Dread(dsetid,memtype,mspace,dspace,H5P_DEFAULT,(void*)&values[0]);
400  if (read_ret < 0) {
401  H5Tclose(memtype);
402  H5Tclose(dtypeid);
403  H5Sclose(dspace);
404  throw InternalErr (__FILE__, __LINE__, "Fail to read the HDF5 compound datatype dataset.");
405  }
406 
407  H5Sclose(dspace);
408  has_values = true;
409  } // end of "if(false == has_values)" block
410 
411  HDF5Structure *h5s = NULL;
412  hid_t memb_id = -1;
413  char* memb_name = NULL;
414 
415  try {
416 
417  // Loop through all the elements in this compound datatype array.
418  for (int element = 0; element < nelms; ++element) {
419 
420  h5s = dynamic_cast<HDF5Structure*>(var()->ptr_duplicate());
421  H5T_class_t memb_cls = H5T_NO_CLASS;
422  int nmembs = 0;
423  size_t memb_offset = 0;
424 
425  if((nmembs = H5Tget_nmembers(memtype)) < 0)
426  throw InternalErr (__FILE__, __LINE__, "Fail to obtain number of HDF5 compound datatype.");
427 
428  for(unsigned int u = 0; u < (unsigned)nmembs; u++) {
429 
430  // Get member type ID
431  if((memb_id = H5Tget_member_type(memtype, u)) < 0)
432  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the datatype of an HDF5 compound datatype member.");
433 
434  // Get member type class
435  if((memb_cls = H5Tget_member_class (memtype, u)) < 0)
436  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the datatype class of an HDF5 compound datatype member.");
437 
438  // Get member offset,H5Tget_member_offset only fails
439  // when H5Tget_memeber_class fails. Sinc H5Tget_member_class
440  // is checked above. So no need to check the return value.
441  memb_offset= H5Tget_member_offset(memtype,u);
442 
443  // Get member name
444  memb_name = H5Tget_member_name(memtype,u);
445  if(memb_name == NULL)
446  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the name of an HDF5 compound datatype member.");
447 
448  BaseType *field = h5s->var(memb_name);
449  if (memb_cls == H5T_COMPOUND) {
450  HDF5Structure &memb_h5s = dynamic_cast<HDF5Structure&>(*field);
451  memb_h5s.do_structure_read(dsetid, memb_id,values,has_values,memb_offset+values_offset+ty_size*element);
452  }
453 
454  else if(memb_cls == H5T_ARRAY) {
455 
456  // memb_id, obtain the number of dimensions
457  int at_ndims = H5Tget_array_ndims(memb_id);
458  if(at_ndims <= 0)
459  throw InternalErr (__FILE__, __LINE__, "Fail to obtain number of dimensions of the array datatype.");
460 
461 #if 0
462  //HDF5Array &h5_array_type = dynamic_cast<HDF5Array&>(*h5s->var(memb_name));
463 #endif
464  HDF5Array &h5_array_type = dynamic_cast<HDF5Array&>(*field);
465  vector<int> at_offset(at_ndims,0);
466  vector<int> at_count(at_ndims,0);
467  vector<int> at_step(at_ndims,0);
468 
469  int at_nelms = h5_array_type.format_constraint(&at_offset[0],&at_step[0],&at_count[0]);
470 
471  // Read the array data
472  h5_array_type.do_h5_array_type_read(dsetid,memb_id,values,has_values,memb_offset+values_offset+ty_size*element,
473  at_nelms,&at_offset[0],&at_count[0],&at_step[0]);
474 
475  }
476  else if(memb_cls == H5T_INTEGER || memb_cls == H5T_FLOAT) {
477 
478  if(true == promote_char_to_short(memb_cls,memb_id)) {
479  void *src = (void*)(&values[0] + (element*ty_size) + values_offset +memb_offset);
480  char val_int8;
481  memcpy(&val_int8,src,1);
482  short val_short=(short)val_int8;
483  field->val2buf(&val_short);
484  }
485  else {
486  field->val2buf(&values[0] + (element*ty_size) + values_offset +memb_offset);
487  }
488 
489  }
490  else if(memb_cls == H5T_STRING) {
491 
492  // distinguish between variable length and fixed length
493  if(true == H5Tis_variable_str(memb_id)) {
494  void *src = (void*)(&values[0]+(element*ty_size)+values_offset + memb_offset);
495  string final_str;
496  get_vlen_str_data((char*)src,final_str);
497  field->val2buf(&final_str);
498  }
499  else {// Obtain fixed-size string value
500  void *src = (void*)(&values[0]+(element*ty_size)+values_offset + memb_offset);
501  vector<char> str_val;
502  size_t memb_size = H5Tget_size(memb_id);
503  if (memb_size == 0) {
504  H5Tclose(memb_id);
505  H5free_memory(memb_name);
506  delete h5s;
507  throw InternalErr (__FILE__, __LINE__,"Fail to obtain the size of HDF5 compound datatype.");
508  }
509  str_val.resize(memb_size);
510  memcpy(&str_val[0],src,memb_size);
511  string temp_string(str_val.begin(),str_val.end());
512  field->val2buf(&temp_string);
513  }
514  }
515  else {
516  H5free_memory(memb_name);
517  H5Tclose(memb_id);
518  delete h5s;
519  throw InternalErr (__FILE__, __LINE__,
520  "Only support the field of compound datatype when the field type class is integer, float, string, array or compound..");
521 
522  }
523 
524  // Close member type ID
525  H5Tclose(memb_id);
526  H5free_memory(memb_name);
527  field->set_read_p(true);
528  } // end "for(unsigned u = 0)"
529  h5s->set_read_p(true);
530  set_vec(element,h5s);
531  delete h5s;
532  } // end "for (int element=0"
533 
534  if(true == has_values) {
535  if(-1 == mspace)
536  throw InternalErr(__FILE__, __LINE__, "memory type and memory space for this compound datatype should be valid.");
537 
538  if(H5Dvlen_reclaim(memtype,mspace,H5P_DEFAULT,(void*)&values[0])<0)
539  throw InternalErr(__FILE__, __LINE__, "Unable to reclaim the compound datatype array.");
540  H5Sclose(mspace);
541 
542  }
543 
544  H5Tclose(dtypeid);
545  H5Tclose(memtype);
546 
547  }
548  catch(...) {
549 
550  if(memb_id != -1)
551  H5Tclose(memb_id);
552  if(memb_name != NULL)
553  H5free_memory(memb_name);
554  if(h5s != NULL)
555  delete h5s;
556  if(true == has_values) {
557  if(H5Dvlen_reclaim(memtype,mspace,H5P_DEFAULT,(void*)(&values[0]))<0) {
558  H5Tclose(memtype);
559  H5Sclose(mspace);
560  }
561  H5Sclose(mspace);
562  }
563  H5Tclose(memtype);
564  H5Tclose(dtypeid);
565  throw;
566  }
567 
568  set_read_p(true);
569 
570  return false;
571 }
572 
573 // Haven't checked the codes and comments
574 // Haven't added the close handles routines for error handlings yet. KY 2011-11-18
575 bool HDF5Array::m_array_of_reference(hid_t dset_id,hid_t dtype_id)
576 {
577 
578 #if (H5_VERS_MAJOR == 1 && (H5_VERS_MINOR == 10 || H5_VERS_MINOR == 8 || H5_VERS_MINOR == 6))
579  hid_t d_dset_id = dset_id;
580  hdset_reg_ref_t *rbuf = NULL;
581 
582  try {
583 
584  vector<int> offset(d_num_dim);
585  vector<int> count(d_num_dim);
586  vector<int> step(d_num_dim);
587 
588 
589  int nelms = format_constraint(&offset[0], &step[0], &count[0]); // Throws Error.
590  vector<string> v_str(nelms);
591 
592  BESDEBUG("h5", "=read() URL type is detected. "
593  << "nelms=" << nelms << " full_size=" << d_num_elm << endl);
594 
595 
596  // Handle regional reference.
597  if (H5Tequal(dtype_id, H5T_STD_REF_DSETREG) < 0) {
598  throw InternalErr(__FILE__, __LINE__, "H5Tequal() failed");
599  }
600 
601  if (H5Tequal(dtype_id, H5T_STD_REF_DSETREG) > 0) {
602  BESDEBUG("h5", "=read() Got regional reference. " << endl);
603  // Vector doesn't work for this case. somehow it doesn't support the type.
604  rbuf = new hdset_reg_ref_t[d_num_elm];
605  if(rbuf == NULL){
606  throw InternalErr(__FILE__, __LINE__, "new() failed.");
607  }
608  if (H5Dread(d_dset_id, H5T_STD_REF_DSETREG, H5S_ALL, H5S_ALL, H5P_DEFAULT, &rbuf[0]) < 0) {
609  throw InternalErr(__FILE__, __LINE__, "H5Dread() failed.");
610  }
611 
612  for (int i = 0; i < nelms; i++) {
613  // Let's assume that URL array is always 1 dimension.
614  BESDEBUG("h5", "=read() rbuf[" << i << "]" <<
615  rbuf[offset[0] + i * step[0]] << endl);
616 
617  if (rbuf[offset[0] + i * step[0]][0] != '\0') {
618  char r_name[DODS_NAMELEN];
619 
620  hid_t did_r = H5RDEREFERENCE(d_dset_id, H5R_DATASET_REGION, rbuf[offset[0] + i * step[0]]);
621  if (did_r < 0) {
622  throw InternalErr(__FILE__, __LINE__, "H5RDEREFERENCE() failed.");
623 
624  }
625 
626  if (H5Iget_name(did_r, (char *) r_name, DODS_NAMELEN) < 0) {
627  throw InternalErr(__FILE__, __LINE__, "H5Iget_name() failed.");
628  }
629  BESDEBUG("h5", "=read() dereferenced name is " << r_name
630  << endl);
631 
632  string varname(r_name);
633  hid_t space_id = H5Rget_region(did_r, H5R_DATASET_REGION, rbuf[offset[0] + i * step[0]]);
634  if (space_id < 0) {
635  throw InternalErr(__FILE__, __LINE__, "H5Rget_region() failed.");
636 
637  }
638 
639  int ndim = H5Sget_simple_extent_ndims(space_id);
640  if (ndim < 0) {
641  throw InternalErr(__FILE__, __LINE__, "H5Sget_simple_extent_ndims() failed.");
642  }
643 
644  BESDEBUG("h5", "=read() dim is " << ndim << endl);
645 
646  string expression;
647  switch (H5Sget_select_type(space_id)) {
648 
649  case H5S_SEL_NONE:
650  BESDEBUG("h5", "=read() None selected." << endl);
651  break;
652 
653  case H5S_SEL_POINTS: {
654  BESDEBUG("h5", "=read() Points selected." << endl);
655  hssize_t npoints = H5Sget_select_npoints(space_id);
656  if (npoints < 0) {
657  throw InternalErr(__FILE__, __LINE__,
658  "Cannot determine number of elements in the dataspace selection");
659  }
660 
661  BESDEBUG("h5", "=read() npoints are " << npoints
662  << endl);
663  vector<hsize_t> buf(npoints * ndim);
664  if (H5Sget_select_elem_pointlist(space_id, 0, npoints, &buf[0]) < 0) {
665  throw InternalErr(__FILE__, __LINE__, "H5Sget_select_elem_pointlist() failed.");
666  }
667 
668 #if 0
669  for (int j = 0; j < npoints * ndim; j++) {
670  "h5", "=read() npoints buf[0] =" << buf[j] <<endl;
671  }
672 #endif
673 
674  for (int j = 0; j < (int) npoints; j++) {
675  // Name of the dataset.
676  expression.append(varname);
677  for (int k = 0; k < ndim; k++) {
678  ostringstream oss;
679  oss << "[" << (int) buf[j * ndim + k] << "]";
680  expression.append(oss.str());
681  }
682  if (j != (int) (npoints - 1)) {
683  expression.append(",");
684  }
685  }
686  v_str[i].append(expression);
687 
688  break;
689  }
690  case H5S_SEL_HYPERSLABS: {
691  vector<hsize_t> start(ndim);
692  vector<hsize_t> end(ndim);
693  vector<hsize_t>stride(ndim);
694  vector<hsize_t>count(ndim);
695  vector<hsize_t>block(ndim);
696 
697  BESDEBUG("h5", "=read() Slabs selected." << endl);
698  BESDEBUG("h5", "=read() nblock is " <<
699  H5Sget_select_hyper_nblocks(space_id) << endl);
700 
701 #if (H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 8)
702  if (H5Sget_select_bounds(space_id, &start[0], &end[0]) < 0) {
703  throw InternalErr(__FILE__, __LINE__, "H5Sget_select_bounds() failed.");
704  }
705 #else
706  if (H5Sget_regular_hyperslab(space_id, &start[0], &stride[0], &count[0], &block[0]) < 0) {
707  throw InternalErr(__FILE__, __LINE__, "H5Sget_regular_hyperslab() failed.");
708  }
709 #endif
710 
711  for (int j = 0; j < ndim; j++) {
712  ostringstream oss;
713  BESDEBUG("h5", "start " << start[j]
714  << "stride "<<stride[j]
715  << "count "<< count[j]
716  << "block "<< block[j]
717  <<endl);
718 
719  // Map from HDF5's start,stride,count,block to DAP's start,stride,end.
720  end[j] = start[j] + stride[j]*(count[j]-1)+(block[j]-1);
721  BESDEBUG("h5", "=read() start is " << start[j]
722  << "=read() end is " << end[j] << endl);
723  oss << "[" << start[j] << ":" << stride[j] << ":" << end[j] << "]";
724  expression.append(oss.str());
725  BESDEBUG("h5", "=read() expression is "
726  << expression << endl)
727  ;
728  }
729  v_str[i] = varname;
730  if (!expression.empty()) {
731  v_str[i].append(expression);
732  }
733  // Constraint expression. [start:1:end]
734  break;
735  }
736  case H5S_SEL_ALL:
737  BESDEBUG("h5", "=read() All selected." << endl);
738  break;
739 
740  default:
741  BESDEBUG("h5", "Unknown space type." << endl);
742  break;
743  }
744 
745  }
746  else {
747  v_str[i] = "";
748  }
749  }
750  delete[] rbuf;
751  }
752 
753  // Handle object reference.
754  if (H5Tequal(dtype_id, H5T_STD_REF_OBJ) < 0) {
755  throw InternalErr(__FILE__, __LINE__, "H5Tequal() failed.");
756  }
757 
758  if (H5Tequal(dtype_id, H5T_STD_REF_OBJ) > 0) {
759  BESDEBUG("h5", "=read() Got object reference. " << endl);
760  vector<hobj_ref_t> orbuf;
761  orbuf.resize(d_num_elm);
762  if (H5Dread(d_dset_id, H5T_STD_REF_OBJ, H5S_ALL, H5S_ALL, H5P_DEFAULT, &orbuf[0]) < 0) {
763  throw InternalErr(__FILE__, __LINE__, "H5Dread failed()");
764  }
765 
766  for (int i = 0; i < nelms; i++) {
767  // Let's assume that URL array is always 1 dimension.
768  hid_t did_r = H5RDEREFERENCE(d_dset_id, H5R_OBJECT, &orbuf[offset[0] + i * step[0]]);
769  if (did_r < 0) {
770  throw InternalErr(__FILE__, __LINE__, "H5RDEREFERENCE() failed.");
771  }
772  char r_name[DODS_NAMELEN];
773  if (H5Iget_name(did_r, (char *) r_name, DODS_NAMELEN) < 0) {
774  throw InternalErr(__FILE__, __LINE__, "H5Iget_name() failed.");
775  }
776 
777  // Shorten the dataset name
778  string varname(r_name);
779 
780  BESDEBUG("h5", "=read() dereferenced name is " << r_name <<endl);
781  v_str[i] = varname;
782  }
783  }
784  set_value(&v_str[0], nelms);
785  return false;
786  }
787  catch (...) {
788  if(rbuf!= NULL)
789  delete[] rbuf;
790  throw;
791  }
792 #else
793  return m_array_of_reference_new_h5_apis(dset_id,dtype_id);
794 
795 #endif
796 }
797 
798 bool HDF5Array::m_array_of_reference_new_h5_apis(hid_t dset_id,hid_t dtype_id) {
799 
800 #if (H5_VERS_MAJOR == 1 && (H5_VERS_MINOR == 10 || H5_VERS_MINOR == 8 || H5_VERS_MINOR == 6))
801  throw InternalErr(__FILE__, __LINE__,
802  "The HDF5 handler compiled with earlier version (<=110)of the HDF5 library should not call method that uses new reference APIs");
803  return false;
804 #else
805 
806  H5R_ref_t *rbuf = NULL;
807  hid_t mem_space_id;
808  hid_t file_space_id;
809 
810  try {
811 
812  // First we need to read the reference data from DAP's hyperslab selection.
813  vector<int> offset(d_num_dim);
814  vector<int> count(d_num_dim);
815  vector<int> step(d_num_dim);
816  vector<hsize_t> hoffset(d_num_dim);
817  vector<hsize_t>hcount(d_num_dim);
818  vector<hsize_t>hstep(d_num_dim);
819 
820  int nelms = format_constraint(&offset[0], &step[0], &count[0]);
821  for (int i = 0; i <d_num_dim; i++) {
822  hoffset[i] = (hsize_t) offset[i];
823  hcount[i] = (hsize_t) count[i];
824  hstep[i] = (hsize_t) step[i];
825  }
826 
827  BESDEBUG("h5", "=read() URL type is detected. "
828  << "nelms=" << nelms << endl);
829 
830  rbuf = new H5R_ref_t[nelms];
831 
832  file_space_id = H5Dget_space(dset_id);
833  if(file_space_id < 0)
834  throw InternalErr(__FILE__, __LINE__, "Fail to obtain reference dataset file space.");
835 
836  if (H5Sselect_hyperslab(file_space_id, H5S_SELECT_SET,
837  &hoffset[0], &hstep[0],
838  &hcount[0], NULL) < 0)
839  throw InternalErr (__FILE__, __LINE__, "Fail to select the hyperslab for reference dataset.");
840 
841 
842  mem_space_id = H5Screate_simple(d_num_dim,&hcount[0],NULL);
843  if(mem_space_id < 0)
844  throw InternalErr(__FILE__, __LINE__, "Fail to obtain reference dataset memory space.");
845 
846  if(H5Dread(dset_id,H5T_STD_REF,mem_space_id,file_space_id,H5P_DEFAULT,&rbuf[0])<0)
847  throw InternalErr(__FILE__, __LINE__, "Fail to read hyperslab reference dataset.");
848 
849  H5Sclose(mem_space_id);
850  H5Sclose(file_space_id);
851 
852  // Now we need to retrieve the reference info. fot the nelms elements.
853  vector<string> v_str;
854 
855  H5R_type_t ref_type = H5Rget_type((const H5R_ref_t *)&rbuf[0]);
856 
857  // The referenced objects can only be either objects or dataset regions.
858  if(ref_type != H5R_OBJECT2 && ref_type !=H5R_DATASET_REGION2)
859  throw InternalErr(__FILE__, __LINE__, "Unsupported reference: neither object nor region references");
860 
861  for (int i = 0; i < nelms; i++) {
862 
863  hid_t obj_id = H5Ropen_object((H5R_ref_t *)&rbuf[i], H5P_DEFAULT, H5P_DEFAULT);
864  if(obj_id < 0)
865  throw InternalErr(__FILE__, __LINE__, "Cannot open the object the reference points to");
866 
867  vector<char> objname;
868  ssize_t objnamelen = -1;
869  if ((objnamelen= H5Iget_name(obj_id,NULL,0))<=0) {
870  H5Oclose(obj_id);
871  throw InternalErr(__FILE__, __LINE__, "Cannot obtain the name length of the object the reference points to");
872  }
873  objname.resize(objnamelen+1);
874  if ((objnamelen= H5Iget_name(obj_id,&objname[0],objnamelen+1))<=0) {
875  H5Oclose(obj_id);
876  throw InternalErr(__FILE__, __LINE__, "Cannot obtain the name length of the object the reference points to");
877  }
878 
879  string objname_str = string(objname.begin(),objname.end());
880  string trim_objname = objname_str.substr(0,objnamelen);
881 
882  // For object references, we just need to save the object full path.
883  if(ref_type == H5R_OBJECT2)
884  v_str.push_back(trim_objname);
885  else {// Must be region reference.
886  H5O_type_t obj_type;
887  if(H5Rget_obj_type3((H5R_ref_t *)&rbuf[i], H5P_DEFAULT, &obj_type) < 0){
888  H5Oclose(obj_id);
889  throw InternalErr(__FILE__, __LINE__, "H5Rget_obj_type3() failed.");
890  }
891  if(obj_type != H5O_TYPE_DATASET) {
892  H5Oclose(obj_id);
893  throw InternalErr(__FILE__, __LINE__, "Region reference must point to a dataset.");
894  }
895  hid_t region_space_id = H5Ropen_region(&rbuf[i],H5P_DEFAULT,H5P_DEFAULT);
896  if (region_space_id < 0) {
897  H5Oclose(obj_id);
898  throw InternalErr(__FILE__, __LINE__, "Cannot obtain the space ID the reference points to");
899  }
900 
901  int ndim = H5Sget_simple_extent_ndims(region_space_id);
902  if (ndim < 0) {
903  H5Sclose(region_space_id);
904  H5Oclose(obj_id);
905  throw InternalErr(__FILE__, __LINE__, "H5Sget_simple_extent_ndims() failed.");
906  }
907 
908  string expression;
909  switch (H5Sget_select_type(region_space_id)) {
910 
911  case H5S_SEL_NONE:
912  BESDEBUG("h5", "=read() None selected." << endl);
913  break;
914 
915  case H5S_SEL_POINTS: {
916  BESDEBUG("h5", "=read() Points selected." << endl);
917  hssize_t npoints = H5Sget_select_npoints(region_space_id);
918  if (npoints < 0) {
919  H5Sclose(region_space_id);
920  H5Oclose(obj_id);
921  throw InternalErr(__FILE__, __LINE__,
922  "Cannot determine number of elements in the dataspace selection");
923  }
924 
925  BESDEBUG("h5", "=read() npoints are " << npoints
926  << endl);
927  vector<hsize_t> buf(npoints * ndim);
928  if (H5Sget_select_elem_pointlist(region_space_id, 0, npoints, &buf[0]) < 0) {
929  H5Sclose(region_space_id);
930  H5Oclose(obj_id);
931  throw InternalErr(__FILE__, __LINE__, "H5Sget_select_elem_pointlist() failed.");
932  }
933 
934 #if 0
935  for (int j = 0; j < npoints * ndim; j++) {
936  "h5", "=read() npoints buf[0] =" << buf[j] <<endl;
937  }
938 #endif
939 
940  for (int j = 0; j < (int) npoints; j++) {
941  // Name of the dataset.
942  expression.append(trim_objname);
943  for (int k = 0; k < ndim; k++) {
944  ostringstream oss;
945  oss << "[" << (int) buf[j * ndim + k] << "]";
946  expression.append(oss.str());
947  }
948  if (j != (int) (npoints - 1)) {
949  expression.append(",");
950  }
951  }
952  v_str.push_back(expression);
953 
954  break;
955  }
956  case H5S_SEL_HYPERSLABS: {
957  vector<hsize_t> start(ndim);
958  vector<hsize_t> end(ndim);
959  vector<hsize_t>stride(ndim);
960  vector<hsize_t>count(ndim);
961  vector<hsize_t>block(ndim);
962 
963  BESDEBUG("h5", "=read() Slabs selected." << endl);
964  BESDEBUG("h5", "=read() nblock is " <<
965  H5Sget_select_hyper_nblocks(region_space_id) << endl);
966 
967  if (H5Sget_regular_hyperslab(region_space_id, &start[0], &stride[0], &count[0], &block[0]) < 0) {
968  H5Sclose(region_space_id);
969  H5Oclose(obj_id);
970  throw InternalErr(__FILE__, __LINE__, "H5Sget_regular_hyperslab() failed.");
971  }
972 
973  expression.append(trim_objname);
974  for (int j = 0; j < ndim; j++) {
975  ostringstream oss;
976  BESDEBUG("h5", "start " << start[j]
977  << "stride "<<stride[j]
978  << "count "<< count[j]
979  << "block "<< block[j]
980  <<endl);
981 
982  // Map from HDF5's start,stride,count,block to DAP's start,stride,end.
983  end[j] = start[j] + stride[j]*(count[j]-1)+(block[j]-1);
984  BESDEBUG("h5", "=read() start is " << start[j]
985  << "=read() end is " << end[j] << endl);
986  oss << "[" << start[j] << ":" << stride[j] << ":" << end[j] << "]";
987  expression.append(oss.str());
988  BESDEBUG("h5", "=read() expression is "
989  << expression << endl)
990  ;
991  }
992  v_str.push_back(expression);
993  // Constraint expression. [start:stride:end]
994  break;
995  }
996  case H5S_SEL_ALL:
997  BESDEBUG("h5", "=read() All selected." << endl);
998  break;
999 
1000  default:
1001  BESDEBUG("h5", "Unknown space type." << endl);
1002  break;
1003  }
1004  H5Sclose(region_space_id);
1005  }
1006  H5Oclose(obj_id);
1007  }
1008  for (int i = 0; i<nelms; i++)
1009  H5Rdestroy(&rbuf[i]);
1010  delete[] rbuf;
1011  H5Sclose(mem_space_id);
1012  H5Sclose(file_space_id);
1013  set_value(&v_str[0], nelms);
1014  return false;
1015  }
1016  catch (...) {
1017  if(rbuf!= NULL)
1018  delete[] rbuf;
1019  H5Sclose(mem_space_id);
1020  H5Sclose(file_space_id);
1021  throw;
1022  }
1023 #endif
1024 }
1025 
1026  // Unused code, keep it for a while.
1027  // Decide not check if the datatype is reference or not. This should be done by URL.
1028 #if 0
1029  bool is_new_ref = false;
1030  bool is_obj_ref = false;
1031  bool is_reg_ref = false;
1032 
1033  htri_t ref_ret = H5Tequal(dtype_id, H5T_STD_REF);
1034 
1035  // Check if this is a new reference(>= HDF5 1.12).
1036  if(ref_ret < 0)
1037  throw InternalErr(__FILE__, __LINE__, "H5Tequal() failed to compare H5T_STD_REF");
1038  else if (ref_ret > 0)
1039  is_new_ref = true;
1040 
1041  if(false == is_new_ref) {
1042 
1043  ref_ret = H5Tequal(dtype_id,H5T_STD_REF_OBJ);
1044  // Check if this is the old object reference(< HDF5 1.12).
1045  if(ref_ret < 0)
1046  throw InternalErr(__FILE__, __LINE__, "H5Tequal() failed to compare H5T_STD_REF_OBJ");
1047  else if (ref_ret > 0)
1048  is_obj_ref = true;
1049 
1050  if(ref_ret == 0) {
1051  // Check if this is the old region reference(< HDF5 1.12)
1052  ref_ret = H5Tequal(dtype_id,H5T_STD_REF_DSETREG);
1053  if(ref_ret < 0)
1054  throw InternalErr(__FILE__, __LINE__, "H5Tequal() failed to compare H5T_STD_REF_DSETREG");
1055  else if (ref_ret > 0)
1056  is_reg_ref = true;
1057  }
1058  }
1059  else {
1060  H5R_ref_t ref_type;
1061  ......
1062  }
1063  if(is_obj_ref == false && is_reg_ref == false)
1064  throw InternalErr(__FILE__, __LINE__, "datatype must be either object ref. or region ref.");
1065 #endif
1066 
1067 
1068 void HDF5Array::m_intern_plain_array_data(char *convbuf,hid_t memtype)
1069 {
1070  if (check_h5str(memtype)) {
1071  vector<string> v_str(d_num_elm);
1072  size_t elesize = H5Tget_size(memtype);
1073  if (elesize == 0) {
1074  throw InternalErr(__FILE__, __LINE__, "H5Tget_size() failed.");
1075  }
1076  vector<char> strbuf(elesize + 1);
1077  BESDEBUG("h5", "=read()<check_h5str() element size=" << elesize
1078  << " d_num_elm=" << d_num_elm << endl);
1079 
1080  for (int strindex = 0; strindex < d_num_elm; strindex++) {
1081  get_strdata(strindex, &convbuf[0], &strbuf[0], elesize);
1082  BESDEBUG("h5", "=read()<get_strdata() strbuf=" << &strbuf[0] << endl);
1083  v_str[strindex] = &strbuf[0];
1084  }
1085  set_read_p(true);
1086  val2buf((void *) &v_str[0]);
1087  }
1088  else {
1089  set_read_p(true);
1090  val2buf((void *) convbuf);
1091  }
1092 }
1093 
1094 
1095 bool HDF5Array::do_h5_array_type_read(hid_t dsetid, hid_t memb_id,vector<char>&values,bool has_values,int values_offset,
1096  int at_nelms,int* at_offset,int* at_count, int* at_step){
1097  //1. Call do array first(datatype must be derived) and the value must be set. We don't support Array datatype
1098  // unless it is inside a compound datatype
1099  if(has_values != true)
1100  throw InternalErr (__FILE__, __LINE__, "Only support the retrieval of HDF5 Array datatype values from the parent compound datatype read.");
1101 
1102  hid_t at_base_type = H5Tget_super(memb_id);
1103  if(at_base_type < 0) {
1104  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the basetype of the array datatype.");
1105  }
1106 
1107  // memb_id, obtain the number of dimensions
1108  int at_ndims = H5Tget_array_ndims(memb_id);
1109  if(at_ndims <= 0) {
1110  H5Tclose(at_base_type);
1111  throw InternalErr (__FILE__, __LINE__, "Fail to obtain number of dimensions of the array datatype.");
1112  }
1113 
1114  vector<hsize_t>at_dims_h(at_ndims,0);
1115 
1116  // Obtain the number of elements for each dims
1117  if(H5Tget_array_dims(memb_id,&at_dims_h[0])<0) {
1118  H5Tclose(at_base_type);
1119  throw InternalErr (__FILE__, __LINE__, "Fail to obtain dimensions of the array datatype.");
1120  }
1121  vector<int>at_dims(at_ndims,0);
1122  for(int i = 0;i<at_ndims;i++) {
1123  at_dims[i] = (int)at_dims_h[i];
1124  }
1125  int at_total_nelms = 1;
1126  for (int i = 0; i <at_ndims; i++)
1127  at_total_nelms = at_total_nelms*at_dims[i];
1128 
1129  H5T_class_t array_cls = H5Tget_class(at_base_type);
1130  if(H5T_NO_CLASS == array_cls) {
1131  H5Tclose(at_base_type);
1132  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the datatype class of the array base type.");
1133  }
1134 
1135  size_t at_base_type_size = H5Tget_size(at_base_type);
1136  if(0 == at_base_type_size){
1137  H5Tclose(at_base_type);
1138  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the size of the array base type.");
1139  }
1140 
1141  // H5 Array type, the basetype is COMPOUND.
1142  if(H5T_COMPOUND == array_cls) {
1143 
1144  // These vectors are used to handle subset of array datatype
1145  vector<int>at_end(at_ndims,0);
1146  vector<int>at_pos(at_ndims,0);
1147  for (int i = 0; i< at_ndims; i++){
1148  at_pos[i] = at_offset[i];
1149  at_end[i] = at_offset[i] + (at_count[i] -1)*at_step[i];
1150  }
1151 
1152  int at_orig_index = INDEX_nD_TO_1D(at_dims,at_pos);
1153 
1154  // To read the array of compound (structure) in DAP, one must read one element each. set_vec is used afterwards.
1155  for (int array_index = 0; array_index <at_nelms; array_index++) {
1156 
1157  // The basetype of the array datatype is compound,-- check if the following line is valid.
1158  HDF5Structure *h5s = dynamic_cast<HDF5Structure*>(var()->ptr_duplicate());
1159  hid_t child_memb_id;
1160  H5T_class_t child_memb_cls;
1161  int child_nmembs;
1162  size_t child_memb_offset;
1163  unsigned child_u;
1164 
1165  if((child_nmembs = H5Tget_nmembers(at_base_type)) < 0) {
1166  H5Tclose(at_base_type);
1167  delete h5s;
1168  throw InternalErr (__FILE__, __LINE__, "Fail to obtain number of HDF5 compound datatype.");
1169  }
1170 
1171  for(child_u = 0; child_u < (unsigned)child_nmembs; child_u++) {
1172 
1173  // Get member type ID
1174  if((child_memb_id = H5Tget_member_type(at_base_type, child_u)) < 0) {
1175  H5Tclose(at_base_type);
1176  delete h5s;
1177  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the datatype of an HDF5 compound datatype member.");
1178  }
1179 
1180  // Get member type class
1181  if((child_memb_cls = H5Tget_member_class (at_base_type, child_u)) < 0) {
1182  H5Tclose(child_memb_id);
1183  H5Tclose(at_base_type);
1184  delete h5s;
1185  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the datatype class of an HDF5 compound datatype member.");
1186  }
1187 
1188  // Get member offset
1189  child_memb_offset= H5Tget_member_offset(at_base_type,child_u);
1190 
1191  // Get member name
1192  char *child_memb_name = H5Tget_member_name(at_base_type,child_u);
1193  if(child_memb_name == NULL) {
1194  H5Tclose(child_memb_id);
1195  H5Tclose(at_base_type);
1196  delete h5s;
1197  throw InternalErr (__FILE__, __LINE__, "Fail to obtain the name of an HDF5 compound datatype member.");
1198  }
1199 
1200  BaseType *field = h5s->var(child_memb_name);
1201  if (child_memb_cls == H5T_COMPOUND) {
1202 
1203  HDF5Structure &memb_h5s = dynamic_cast<HDF5Structure&>(*field);
1204  //
1205  // Call structure read when reading the whole array. sa1{sa2[100]}
1206  // sa2[100] is an array datatype.
1207  // If reading the whole buffer, just add the buffer.
1208  if(at_total_nelms == at_nelms) {
1209  memb_h5s.do_structure_read(dsetid,child_memb_id, values,has_values,values_offset+at_base_type_size*array_index+child_memb_offset);
1210  }
1211  // Subset of sa2, sa2[10:100:2]; sa2[100] is an array datatype. The whole array sa2[100] is to be read into somewhere in buffer values.
1212  else {// The subset should be considered. adjust memb_offset+values_offset+???,make sure only the subset is selected.
1213  // at_total_nelms is 100 but at_nelms is (100-10)/2+1=46. The starting point of the whole array is values+memb_offset_values_offset
1214  // When the datatype is structure, we have to obtain the index one by one.
1215 
1216  memb_h5s.do_structure_read(dsetid, child_memb_id, values,has_values,values_offset+at_base_type_size*at_orig_index+child_memb_offset);
1217 
1218  }
1219 
1220  }
1221  else if(child_memb_cls == H5T_ARRAY) {
1222 
1223  // memb_id, obtain the number of dimensions
1224  int child_at_ndims = H5Tget_array_ndims(child_memb_id);
1225  if(child_at_ndims <= 0) {
1226  H5Tclose(at_base_type);
1227  H5Tclose(child_memb_id);
1228  H5free_memory(child_memb_name);
1229  delete h5s;
1230  throw InternalErr (__FILE__, __LINE__, "Fail to obtain number of dimensions of the array datatype.");
1231 
1232  }
1233 
1234  HDF5Array &h5_array_type = dynamic_cast<HDF5Array&>(*field);
1235  vector<int> child_at_offset(child_at_ndims,0);
1236  vector<int> child_at_count(child_at_ndims,0);
1237  vector<int> child_at_step(child_at_ndims,0);
1238 
1239  int child_at_nelms = h5_array_type.format_constraint(&child_at_offset[0],&child_at_step[0],&child_at_count[0]);
1240 
1241  if(at_total_nelms == at_nelms) {
1242  h5_array_type.do_h5_array_type_read(dsetid,child_memb_id,values,has_values,child_memb_offset+values_offset+at_base_type_size*array_index,
1243  child_at_nelms,&child_at_offset[0],&child_at_count[0],&child_at_step[0]);
1244  }
1245  else {// Adjust memb_offset+values_offset, basically change at_base_type_size*array_index
1246  h5_array_type.do_h5_array_type_read(dsetid,child_memb_id,values,has_values,child_memb_offset+values_offset+at_base_type_size*at_orig_index,
1247  child_at_nelms,&child_at_offset[0],&child_at_count[0],&child_at_step[0]);
1248 
1249  }
1250  }
1251 
1252  else if(H5T_INTEGER == child_memb_cls || H5T_FLOAT == child_memb_cls){
1253 
1254  int number_index =((at_total_nelms == at_nelms)?array_index:at_orig_index);
1255  if(true == promote_char_to_short(child_memb_cls,child_memb_id)) {
1256  void *src = (void*)(&values[0] + (number_index*at_base_type_size) + values_offset +child_memb_offset);
1257  char val_int8;
1258  memcpy(&val_int8,src,1);
1259  short val_short=(short)val_int8;
1260  field->val2buf(&val_short);
1261  }
1262  else
1263  field->val2buf(&values[0] + (number_index * at_base_type_size) + values_offset+child_memb_offset);
1264 
1265 
1266  }
1267  else if(H5T_STRING == child_memb_cls){
1268 
1269  int string_index =((at_total_nelms == at_nelms)?array_index:at_orig_index);
1270 
1271  // distinguish between variable length and fixed length
1272  if(true == H5Tis_variable_str(child_memb_id)) {
1273 
1274  // Need to check if the size of variable length array type is right in HDF5 lib.
1275  void *src = (void*)(&values[0]+(string_index *at_base_type_size)+values_offset+child_memb_offset);
1276  string final_str;
1277  char*temp_bp =(char*)src;
1278  get_vlen_str_data(temp_bp,final_str);
1279  field->val2buf(&final_str[0]); //field->set_value(final_str);
1280 
1281  }
1282  else {// Obtain string
1283  void *src = (void*)(&values[0]+(string_index *at_base_type_size)+values_offset+child_memb_offset);
1284  vector<char> str_val;
1285  size_t memb_size = H5Tget_size(child_memb_id);
1286  if (memb_size == 0) {
1287  H5Tclose(child_memb_id);
1288  H5Tclose(at_base_type);
1289  H5free_memory(child_memb_name);
1290  delete h5s;
1291  throw InternalErr (__FILE__, __LINE__,"Fail to obtain the size of HDF5 compound datatype.");
1292  }
1293  str_val.resize(memb_size);
1294  memcpy(&str_val[0],src,memb_size);
1295  field->val2buf(&str_val[0]);
1296 
1297  }
1298  }
1299  else {
1300  H5Tclose(child_memb_id);
1301  H5Tclose(at_base_type);
1302  H5free_memory(child_memb_name);
1303  delete h5s;
1304  throw InternalErr (__FILE__, __LINE__, "Unsupported datatype class for the array base type.");
1305 
1306 
1307  }
1308  field->set_read_p(true);
1309  H5free_memory(child_memb_name);
1310  H5Tclose(child_memb_id);
1311 
1312  } // end "for ( child_u = 0)"
1313  h5s->set_read_p(true);
1314 
1315  // Save the value of this element to DAP structure.
1316  set_vec(array_index,h5s);
1317  delete h5s;
1318 
1319  vector<int>at_offsetv(at_pos.size(),0);
1320  vector<int>at_stepv(at_pos.size(),0);
1321  for (unsigned int at_index = 0; at_index<at_pos.size();at_index++){
1322  at_offsetv[at_index] = at_offset[at_index];
1323  at_stepv[at_index] = at_step[at_index];
1324  }
1325  //obtain the next position of the selected point based on the offset,end and step.
1326  obtain_next_pos(at_pos,at_offsetv,at_end,at_stepv,(int)(at_pos.size()));
1327  at_orig_index = INDEX_nD_TO_1D(at_dims,at_pos);
1328  }// end for "(array_index = 0) for array (compound)datatype"
1329 
1330  // Need to check the usage of set_read_p(true);
1331 #if 0
1332  //set_read_p(true);
1333 #endif
1334  }
1335  else if(H5T_INTEGER == array_cls|| H5T_FLOAT == array_cls) {
1336 
1337  // If no subset for the array datatype, just read the whole buffer.
1338  if(at_total_nelms == at_nelms) {
1339 
1340  // For DAP2 char should be mapped to short
1341  if( true == promote_char_to_short(array_cls ,at_base_type)) {
1342  vector<char> val_int8;
1343  val_int8.resize(at_nelms);
1344  void*src = (void*)(&values[0] +values_offset);
1345  memcpy(&val_int8[0],src,at_nelms);
1346 
1347  vector<short> val_short;
1348  for (int i = 0; i<at_nelms; i++)
1349  val_short[i] = (short)val_int8[i];
1350 
1351  val2buf(&val_short[0]);
1352 
1353  }
1354  else // short cut for others
1355  val2buf(&values[0] + values_offset);
1356 
1357  }
1358  else { // Adjust the value for the subset of the array datatype
1359 
1360  // Obtain the correponding DAP type of the HDF5 data type
1361  string dap_type = get_dap_type(at_base_type,is_dap4());
1362 
1363  // The total array type data is read.
1364  void*src = (void*)(&values[0] + values_offset);
1365 
1366  // set the original position to the starting point
1367  vector<int>at_pos(at_ndims,0);
1368  for (int i = 0; i< at_ndims; i++)
1369  at_pos[i] = at_offset[i];
1370 
1371  if( BYTE == dap_type) {
1372 
1373  vector<unsigned char>total_val;
1374  total_val.resize(at_total_nelms);
1375  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1376 
1377  vector<unsigned char>final_val;
1378  subset<unsigned char>(
1379  &total_val[0],
1380  at_ndims,
1381  at_dims,
1382  at_offset,
1383  at_step,
1384  at_count,
1385  &final_val,
1386  at_pos,
1387  0
1388  );
1389 
1390  set_value((dods_byte*)&final_val[0],at_nelms);
1391 
1392 
1393  }
1394  else if( INT16 == dap_type) {
1395 
1396  // promote char to short,DAP2 doesn't have "char" type
1397  if(true == promote_char_to_short(array_cls,at_base_type)) {
1398  vector<char>total_val;
1399  total_val.resize(at_total_nelms);
1400  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1401 
1402  vector<char>final_val;
1403  subset<char>(
1404  &total_val[0],
1405  at_ndims,
1406  at_dims,
1407  at_offset,
1408  at_step,
1409  at_count,
1410  &final_val,
1411  at_pos,
1412  0
1413  );
1414 
1415  vector<short> final_val_short;
1416  final_val_short.resize(at_nelms);
1417  for(int i = 0; i<at_nelms; i++)
1418  final_val_short[i] = final_val[i];
1419 
1420  val2buf(&final_val_short[0]);
1421 
1422  }
1423  else {// short
1424 
1425  vector<short>total_val;
1426  total_val.resize(at_total_nelms);
1427  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1428 
1429  vector<short>final_val;
1430  subset<short>(
1431  &total_val[0],
1432  at_ndims,
1433  at_dims,
1434  at_offset,
1435  at_step,
1436  at_count,
1437  &final_val,
1438  at_pos,
1439  0
1440  );
1441 
1442  val2buf(&final_val[0]);
1443 
1444  }
1445  }
1446  else if( UINT16 == dap_type) {
1447  vector<unsigned short>total_val;
1448  total_val.resize(at_total_nelms);
1449  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1450 
1451  vector<unsigned short>final_val;
1452  subset<unsigned short>(
1453  &total_val[0],
1454  at_ndims,
1455  at_dims,
1456  at_offset,
1457  at_step,
1458  at_count,
1459  &final_val,
1460  at_pos,
1461  0
1462  );
1463 
1464  val2buf(&final_val[0]);
1465 
1466  }
1467  else if(UINT32 == dap_type) {
1468  vector<unsigned int>total_val;
1469  total_val.resize(at_total_nelms);
1470  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1471 
1472  vector<unsigned int>final_val;
1473  subset<unsigned int>(
1474  &total_val[0],
1475  at_ndims,
1476  at_dims,
1477  at_offset,
1478  at_step,
1479  at_count,
1480  &final_val,
1481  at_pos,
1482  0
1483  );
1484  val2buf(&final_val[0]);
1485 
1486 
1487  }
1488  else if(INT32 == dap_type) {
1489  vector<int>total_val;
1490  total_val.resize(at_total_nelms);
1491  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1492 
1493  vector<int>final_val;
1494  subset<int>(
1495  &total_val[0],
1496  at_ndims,
1497  at_dims,
1498  at_offset,
1499  at_step,
1500  at_count,
1501  &final_val,
1502  at_pos,
1503  0
1504  );
1505 
1506  val2buf(&final_val[0]);
1507 
1508  }
1509  else if(FLOAT32 == dap_type) {
1510  vector<float>total_val;
1511  total_val.resize(at_total_nelms);
1512  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1513 
1514  vector<float>final_val;
1515  subset<float>(
1516  &total_val[0],
1517  at_ndims,
1518  at_dims,
1519  at_offset,
1520  at_step,
1521  at_count,
1522  &final_val,
1523  at_pos,
1524  0
1525  );
1526 
1527  val2buf(&final_val[0]);
1528 
1529  }
1530  else if(FLOAT64 == dap_type) {
1531  vector<double>total_val;
1532  total_val.resize(at_total_nelms);
1533  memcpy(&total_val[0],src,at_total_nelms*at_base_type_size);
1534 
1535  vector<double>final_val;
1536  subset<double>(
1537  &total_val[0],
1538  at_ndims,
1539  at_dims,
1540  at_offset,
1541  at_step,
1542  at_count,
1543  &final_val,
1544  at_pos,
1545  0
1546  );
1547 #if 0
1548  //field->val2buf(&final_val[0]);
1549 #endif
1550  val2buf(&final_val[0]);
1551 
1552  }
1553  else {
1554  H5Tclose(at_base_type);
1555  throw InternalErr (__FILE__, __LINE__,
1556  "Non-supported integer or float datatypes");
1557  }
1558 
1559  }
1560  }
1561  else if(H5T_STRING == array_cls) {
1562 
1563  // set the original position to the starting point
1564  vector<int>at_pos(at_ndims,0);
1565  for (int i = 0; i< at_ndims; i++)
1566  at_pos[i] = at_offset[i];
1567 
1568  vector<string>total_strval;
1569  total_strval.resize(at_total_nelms);
1570 
1571  if(true == H5Tis_variable_str(at_base_type)) {
1572  void *src = (void*)(&values[0]+values_offset);
1573  char*temp_bp =(char*)src;
1574  for(int i = 0;i <at_total_nelms; i++){
1575  string tempstrval;
1576  get_vlen_str_data(temp_bp,tempstrval);
1577  total_strval[i] = tempstrval;
1578  temp_bp += at_base_type_size;
1579  }
1580  if(at_total_nelms == at_nelms) {
1581 #if 0
1582  //field->set_value(total_strval,at_total_nelms);
1583 #endif
1584  set_value(total_strval,at_total_nelms);
1585  }
1586  else {// obtain subset for variable-length string.
1587 #if 0
1588  //field->val2buf(&values[0] + values_offset);
1589 #endif
1590  vector<string>final_val;
1591  subset<string>(
1592  &total_strval[0],
1593  at_ndims,
1594  at_dims,
1595  at_offset,
1596  at_step,
1597  at_count,
1598  &final_val,
1599  at_pos,
1600  0
1601  );
1602 
1603  set_value(final_val,at_nelms);
1604 
1605  }
1606 
1607  }
1608  else {// For fixed-size string.
1609  void *src = (void*)(&values[0]+values_offset);
1610  for(int i = 0; i <at_total_nelms; i++)
1611  total_strval[i].resize(at_base_type_size);
1612 
1613  vector<char> str_val;
1614  str_val.resize(at_total_nelms*at_base_type_size);
1615  memcpy((void*)&str_val[0],src,at_total_nelms*at_base_type_size);
1616  string total_in_one_string(str_val.begin(),str_val.end());
1617  for(int i = 0; i<at_total_nelms;i++)
1618  total_strval[i] = total_in_one_string.substr(i*at_base_type_size,at_base_type_size);
1619 
1620  if(at_total_nelms == at_nelms)
1621  set_value(total_strval,at_total_nelms);
1622  else {
1623  vector<string>final_val;
1624  subset<string>(
1625  &total_strval[0],
1626  at_ndims,
1627  at_dims,
1628  at_offset,
1629  at_step,
1630  at_count,
1631  &final_val,
1632  at_pos,
1633  0
1634  );
1635  set_value(final_val,at_nelms);
1636 
1637  }
1638  }
1639 
1640  }
1641  else {
1642  H5Tclose(at_base_type);
1643  throw InternalErr (__FILE__, __LINE__,
1644  "Only support the field of compound datatype when the field type class is integer, float, string, array or compound..");
1645 
1646  }
1647 
1648  H5Tclose(at_base_type);
1649 
1650  return true;
1651 }
1652 
1654 inline int
1655 HDF5Array::INDEX_nD_TO_1D (const std::vector < int > &dims,
1656  const std::vector < int > &pos)
1657 {
1658  //
1659  // "int a[10][20][30] // & a[1][2][3] == a + (20*30+1 + 30*2 + 1 *3)"
1660  // "int b[10][2]; // &b[1][2] == b + (20*1 + 2);"
1661  //
1662  assert (dims.size () == pos.size ());
1663  int sum = 0;
1664  int start = 1;
1665 
1666  for (unsigned int p = 0; p < pos.size (); p++) {
1667  int m = 1;
1668 
1669  for (unsigned int j = start; j < dims.size (); j++)
1670  m *= dims[j];
1671  sum += m * pos[p];
1672  start++;
1673  }
1674  return sum;
1675 }
1676 
1677 // Obtain the dimension index of the next pos. of the point based on the offset, step and end
1678 bool HDF5Array::obtain_next_pos(vector<int>& pos, vector<int>&start,vector<int>&end,vector<int>&step,int rank_change) {
1679 
1680  if((pos[rank_change-1] + step[rank_change-1])<=end[rank_change-1]) {
1681  pos[rank_change-1] = pos[rank_change-1] + step[rank_change-1];
1682  return true;
1683  }
1684  else {
1685  if( 1 == rank_change)
1686  return false;
1687  pos[rank_change-1] = start[rank_change-1];
1688  obtain_next_pos(pos,start,end,step,rank_change-1);
1689  }
1690  return true;
1691 }
1692 
1694 //
1695 // \param input Input variable
1696 // \param dim dimension info of the input
1697 // \param start start indexes of each dim
1698 // \param stride stride of each dim
1699 // \param edge count of each dim
1700 // \param poutput output variable
1701 // \parrm index dimension index
1702 // \return 0 if successful. -1 otherwise.
1703 //
1704 template<typename T>
1705 int HDF5Array::subset(
1706  const T input[],
1707  int rank,
1708  vector<int> & dim,
1709  int start[],
1710  int stride[],
1711  int edge[],
1712  std::vector<T> *poutput,
1713  vector<int>& pos,
1714  int index)
1715 {
1716  for(int k=0; k<edge[index]; k++)
1717  {
1718  pos[index] = start[index] + k*stride[index];
1719  if(index+1<rank)
1720  subset(input, rank, dim, start, stride, edge, poutput,pos,index+1);
1721  if(index==rank-1)
1722  {
1723  poutput->push_back(input[INDEX_nD_TO_1D( dim, pos)]);
1724  }
1725  } // end of for
1726  return 0;
1727 } // end of template<typename T> static int subset
1728 
1729 
1730 // public functions to set all parameters needed in read function.
1731 
1732 
1733 void HDF5Array::set_memneed(size_t need) {
1734  d_memneed = need;
1735 }
1736 
1737 void HDF5Array::set_numdim(int ndims) {
1738  d_num_dim = ndims;
1739 }
1740 
1741 void HDF5Array::set_numelm(int nelms) {
1742  d_num_elm = nelms;
1743 }
1744 
1745 hid_t HDF5Array::mkstr(int size, H5T_str_t pad)
1746 {
1747 
1748  hid_t str_type;
1749 
1750  if ((str_type = H5Tcopy(H5T_C_S1)) < 0)
1751  return -1;
1752  if (H5Tset_size(str_type, (size_t) size) < 0)
1753  return -1;
1754  if (H5Tset_strpad(str_type, pad) < 0)
1755  return -1;
1756 
1757  return str_type;
1758 }
1759 
1760 // We don't inherit libdap Array Class's transform_to_dap4 method since CF option is still using it.
1761 BaseType* HDF5Array::h5dims_transform_to_dap4(D4Group *grp,const vector<string> &dimpath) {
1762 
1763  BESDEBUG("h5", "<h5dims_transform_to_dap4" << endl);
1764 
1765  if(grp == NULL)
1766  return NULL;
1767 
1768  Array *dest = static_cast<HDF5Array*>(ptr_duplicate());
1769 
1770  // If there is just a size, don't make
1771  // a D4Dimension (In DAP4 you cannot share a dimension unless it has
1772  // a name). jhrg 3/18/14
1773 
1774  int k = 0;
1775  for (Array::Dim_iter d = dest->dim_begin(), e = dest->dim_end(); d != e; ++d) {
1776 
1777  if (false == (*d).name.empty()) {
1778  BESDEBUG("h5", "<coming to the dimension loop, has name " << (*d).name<<endl);
1779  BESDEBUG("h5", "<coming to the dimension loop, has dimpath " << dimpath[k] <<endl);
1780  BESDEBUG("h5", "<coming to the dimension loop, has dimpath group " << dimpath[k].substr(0,dimpath[k].find_last_of("/")+1) <<endl);
1781 
1782  D4Group *temp_grp = grp;
1783  D4Dimension *d4_dim = NULL;
1784  bool is_dim_nonc4_grp = false;
1785 
1786  while(temp_grp) {
1787 
1788  BESDEBUG("h5", "<coming to the group has name " << temp_grp->name()<<endl);
1789  BESDEBUG("h5", "<coming to the group has fullpath " << temp_grp->FQN()<<endl);
1790 
1791  //Obtain all the dimensions of this group.
1792  D4Dimensions *temp_dims = temp_grp->dims();
1793 
1794  // Check if this dimension is defined in this group
1795  d4_dim = temp_dims->find_dim((*d).name);
1796 
1797  // Need the full path of the dimension name
1798  string d4_dim_path = dimpath[k].substr(0,dimpath[k].find_last_of("/")+1);
1799  BESDEBUG("h5", "d4_dim_path is " << d4_dim_path<<endl);
1800 
1801  bool ancestor_grp = false;
1802 
1803  // If the dim_path is within this group or its ancestor, this is valid.
1804  if(d4_dim_path.find(temp_grp->FQN())==0 || temp_grp->FQN().find(d4_dim_path)==0)
1805  ancestor_grp = true;
1806 
1807  // If we find this dimension and the dimension is on the ancestral path,
1808  // this follows the netCDF-4/DAP4 dimension model, break.
1809  if(d4_dim && (temp_grp->FQN() == d4_dim_path)) {
1810  BESDEBUG("h5", "<FInd dimension name " << (*d).name<<endl);
1811  (*d).dim = d4_dim;
1812  is_dim_nonc4_grp = false;
1813  break;
1814  }
1815  // If the dimension name is not on the ancestral path, this
1816  // dimension must be on another path, mark it.
1817  //else if( ancestor_grp == false && is_dim_nonc4_grp == false) {
1818  else if( ancestor_grp == false) {
1819  is_dim_nonc4_grp = true;
1820  break;
1821  }
1822  else
1823  d4_dim = NULL;
1824 
1825  if(temp_grp->get_parent())
1826  temp_grp = static_cast<D4Group*>(temp_grp->get_parent());
1827  else
1828  temp_grp = 0;
1829 
1830  }
1831 
1832  // Not find this dimension in any of the ancestor groups, add it to this group.
1833  // The following block is fine, but to avoid the complaint from sonarcloud.
1834  // Use a bool.
1835  if(true == is_dim_nonc4_grp) {
1836  string err= "The variable " + var_path +" has dimension ";
1837  err += dimpath[k] + ". This dimension is not under its ancestor or the current group.";
1838  err += " This is not supported.";
1839  delete dest;
1840  throw InternalErr(__FILE__,__LINE__,err);
1841  }
1842 
1843  bool d4_dim_null = ((d4_dim==NULL)?true:false);
1844  if(d4_dim_null == true) {
1845  d4_dim = new D4Dimension((*d).name, (*d).size);
1846  D4Dimensions * dims = grp->dims();
1847  BESDEBUG("h5", "<Just before adding D4 dimension to group" << endl);
1848  dims->add_dim_nocopy(d4_dim);
1849  (*d).dim = d4_dim;
1850  }
1851  }
1852  k++;
1853  }
1854 
1855  dest->set_is_dap4(true);
1856 
1857  return dest;
1858 
1859 }
A class for handling all types of array in HDF5 for the default option.
This class that translates HDF5 string into DAP string for the default option.
This class converts HDF5 compound type into DAP structure for the default option.
virtual libdap::BaseType * ptr_duplicate()
Definition: HDF5Array.cc:54
virtual bool read()
Reads HDF5 array data into local buffer.
Definition: HDF5Array.cc:115
void set_numdim(int ndims)
remembers number of dimensions of this array.
Definition: HDF5Array.cc:1737
HDF5Array(const std::string &n, const std::string &d, libdap::BaseType *v)
Constructor.
Definition: HDF5Array.cc:58
void set_numelm(int nelms)
remembers number of elements in this array.
Definition: HDF5Array.cc:1741
void set_memneed(size_t need)
remembers memory size needed.
Definition: HDF5Array.cc:1733
int get_slabdata(hid_t dset, int *offset, int *step, int *count, int num_dim, void *buf)
Definition: h5common.cc:144
void get_data(hid_t dset, void *buf)
Definition: h5common.cc:50
void get_strdata(int strindex, char *allbuf, char *buf, int elesize)
Definition: h5common.cc:115
bool check_h5str(hid_t h5type)
Definition: h5get.cc:842
string get_dap_type(hid_t type, bool is_dap4)
Definition: h5get.cc:291
const int DODS_NAMELEN
Maximum length of variable or attribute name(default option only).
Definition: hdf5_handler.h:65