NetCDF  4.3.2
nc4file.c
Go to the documentation of this file.
00001 
00011 #include "config.h"
00012 #include <errno.h>  /* netcdf functions sometimes return system errors */
00013 
00014 #include "nc.h"
00015 #include "nc4internal.h"
00016 #include "nc4dispatch.h"
00017 
00018 #include "H5DSpublic.h"
00019 
00020 #ifdef USE_HDF4
00021 #include <mfhdf.h>
00022 #endif
00023 
00024 #if 0 /*def USE_PNETCDF*/
00025 #include <pnetcdf.h>
00026 #endif
00027 
00028 /* This is to track opened HDF5 objects to make sure they are
00029  * closed. */
00030 #ifdef EXTRA_TESTS
00031 extern int num_plists;
00032 extern int num_spaces;
00033 #endif /* EXTRA_TESTS */
00034 
00035 #define MIN_DEFLATE_LEVEL 0
00036 #define MAX_DEFLATE_LEVEL 9
00037 
00038 /* These are the special attributes added by the HDF5 dimension scale
00039  * API. They will be ignored by netCDF-4. */
00040 #define REFERENCE_LIST "REFERENCE_LIST"
00041 #define CLASS "CLASS"
00042 #define DIMENSION_LIST "DIMENSION_LIST"
00043 #define NAME "NAME"
00044 
00045 /* Struct to track information about objects in a group, for nc4_rec_read_metadata() */
00046 typedef struct NC4_rec_read_metadata_obj_info
00047 {
00048     hid_t oid;                          /* HDF5 object ID */
00049     char oname[NC_MAX_NAME + 1];        /* Name of object */
00050     H5G_stat_t statbuf;                 /* Information about the object */
00051     struct NC4_rec_read_metadata_obj_info *next;        /* Pointer to next node in list */
00052 } NC4_rec_read_metadata_obj_info_t;
00053 
00054 /* User data struct for call to H5Literate() in nc4_rec_read_metadata() */
00055 /* Tracks the groups, named datatypes and datasets in the group, for later use */
00056 typedef struct NC4_rec_read_metadata_ud
00057 {
00058    NC4_rec_read_metadata_obj_info_t *grps_head, *grps_tail;      /* Pointers to head & tail of list of groups */
00059    NC_GRP_INFO_T *grp;                                          /* Pointer to parent group */
00060 } NC4_rec_read_metadata_ud_t;
00061 
00062 /* Forward */
00063 static int NC4_enddef(int ncid);
00064 static int nc4_rec_read_metadata(NC_GRP_INFO_T *grp);
00065 static int close_netcdf4_file(NC_HDF5_FILE_INFO_T *h5, int abort);
00066 
00067 /* These are the default chunk cache sizes for HDF5 files created or
00068  * opened with netCDF-4. */
00069 size_t nc4_chunk_cache_size = CHUNK_CACHE_SIZE;
00070 size_t nc4_chunk_cache_nelems = CHUNK_CACHE_NELEMS;
00071 float nc4_chunk_cache_preemption = CHUNK_CACHE_PREEMPTION;
00072 
00073 /* To turn off HDF5 error messages, I have to catch an early
00074    invocation of a netcdf function. */
00075 static int virgin = 1;
00076 
00077 /* For performance, fill this array only the first time, and keep it
00078  * in global memory for each further use. */
00079 #define NUM_TYPES 12
00080 static hid_t h5_native_type_constant_g[NUM_TYPES];
00081 static const char nc_type_name_g[NUM_TYPES][NC_MAX_NAME + 1] = {"char", "byte", "short",
00082     "int", "float", "double", "ubyte",
00083     "ushort", "uint", "int64",
00084     "uint64", "string"};
00085 static const nc_type nc_type_constant_g[NUM_TYPES] = {NC_CHAR, NC_BYTE, NC_SHORT,
00086     NC_INT, NC_FLOAT, NC_DOUBLE, NC_UBYTE,
00087     NC_USHORT, NC_UINT, NC_INT64,
00088     NC_UINT64, NC_STRING};
00089 static const int nc_type_size_g[NUM_TYPES] = {sizeof(char), sizeof(char), sizeof(short), 
00090     sizeof(int), sizeof(float), sizeof(double), sizeof(unsigned char),
00091     sizeof(unsigned short), sizeof(unsigned int), sizeof(long long), 
00092     sizeof(unsigned long long), sizeof(char *)};
00093 
00094 /* Set chunk cache size. Only affects files opened/created *after* it
00095  * is called.  */
00096 int
00097 nc_set_chunk_cache(size_t size, size_t nelems, float preemption)
00098 {
00099    if (preemption < 0 || preemption > 1)
00100       return NC_EINVAL;
00101    nc4_chunk_cache_size = size;
00102    nc4_chunk_cache_nelems = nelems;
00103    nc4_chunk_cache_preemption = preemption;
00104    return NC_NOERR;
00105 }
00106 
00107 /* Get chunk cache size. Only affects files opened/created *after* it
00108  * is called.  */
00109 int
00110 nc_get_chunk_cache(size_t *sizep, size_t *nelemsp, float *preemptionp)
00111 {
00112    if (sizep)
00113       *sizep = nc4_chunk_cache_size;
00114 
00115    if (nelemsp)
00116       *nelemsp = nc4_chunk_cache_nelems;
00117    
00118    if (preemptionp)
00119       *preemptionp = nc4_chunk_cache_preemption;
00120    return NC_NOERR;
00121 }
00122 
00123 /* Required for fortran to avoid size_t issues. */
00124 int
00125 nc_set_chunk_cache_ints(int size, int nelems, int preemption)
00126 {
00127    if (size <= 0 || nelems <= 0 || preemption < 0 || preemption > 100)
00128       return NC_EINVAL;
00129    nc4_chunk_cache_size = size;
00130    nc4_chunk_cache_nelems = nelems;
00131    nc4_chunk_cache_preemption = (float)preemption / 100;
00132    return NC_NOERR;
00133 }
00134 
00135 int
00136 nc_get_chunk_cache_ints(int *sizep, int *nelemsp, int *preemptionp)
00137 {
00138    if (sizep)
00139       *sizep = (int)nc4_chunk_cache_size;
00140    if (nelemsp)
00141       *nelemsp = (int)nc4_chunk_cache_nelems;
00142    if (preemptionp)
00143       *preemptionp = (int)(nc4_chunk_cache_preemption * 100);
00144 
00145    return NC_NOERR;
00146 }
00147 
00148 /* This will return the length of a netcdf data type in bytes. */
00149 int
00150 nc4typelen(nc_type type)
00151 {
00152    switch(type){
00153       case NC_BYTE:
00154       case NC_CHAR:
00155       case NC_UBYTE:
00156          return 1;
00157       case NC_USHORT:
00158       case NC_SHORT:
00159          return 2;
00160       case NC_FLOAT:
00161       case NC_INT:
00162       case NC_UINT:
00163          return 4;
00164       case NC_DOUBLE: 
00165       case NC_INT64:
00166       case NC_UINT64:
00167          return 8;
00168    }
00169    return -1;
00170 }
00171 
00172 /* Given a filename, check to see if it is a HDF5 file. */
00173 #define MAGIC_NUMBER_LEN 4
00174 #define NC_HDF5_FILE 1
00175 #define NC_HDF4_FILE 2
00176 static int
00177 nc_check_for_hdf(const char *path, int use_parallel, MPI_Comm comm, MPI_Info info, 
00178                  int *hdf_file)
00179 {
00180    char blob[MAGIC_NUMBER_LEN];
00181    
00182    assert(hdf_file && path);
00183    LOG((3, "%s: path %s", __func__, path));
00184 
00185    /* HDF5 function handles possible user block at beginning of file */
00186    if(H5Fis_hdf5(path)) 
00187    {
00188        *hdf_file = NC_HDF5_FILE;
00189    } else {
00190 
00191 /* Get the 4-byte blob from the beginning of the file. Don't use posix
00192  * for parallel, use the MPI functions instead. */
00193 #ifdef USE_PARALLEL
00194        if (use_parallel)
00195        {
00196            MPI_File fh;
00197            MPI_Status status;
00198            int retval;
00199            if ((retval = MPI_File_open(comm, (char *)path, MPI_MODE_RDONLY,
00200                                        info, &fh)) != MPI_SUCCESS)
00201                return NC_EPARINIT;
00202            if ((retval = MPI_File_read(fh, blob, MAGIC_NUMBER_LEN, MPI_CHAR,
00203                                        &status)) != MPI_SUCCESS)
00204                return NC_EPARINIT;
00205            if ((retval = MPI_File_close(&fh)) != MPI_SUCCESS)
00206                return NC_EPARINIT;
00207        }
00208        else
00209 #endif /* USE_PARALLEL */
00210        {
00211            FILE *fp;
00212            if (!(fp = fopen(path, "r")) ||
00213                fread(blob, MAGIC_NUMBER_LEN, 1, fp) != 1) {
00214 
00215              if(fp) fclose(fp);
00216              return errno;
00217            }
00218            fclose(fp);
00219        }
00220        
00221        /* Check for HDF4. */
00222        if (!strncmp(blob, "\016\003\023\001", MAGIC_NUMBER_LEN))
00223            *hdf_file = NC_HDF4_FILE;
00224        else
00225            *hdf_file = 0;
00226    }
00227    return NC_NOERR;
00228 }
00229    
00230 /* Create a HDF5/netcdf-4 file. */
00231 
00232 static int
00233 nc4_create_file(const char *path, int cmode, MPI_Comm comm, MPI_Info info,
00234                 NC *nc) 
00235 {
00236    hid_t fcpl_id, fapl_id = -1;
00237    unsigned flags;
00238    FILE *fp;
00239    int retval = NC_NOERR;
00240    NC_HDF5_FILE_INFO_T* nc4_info = NULL;
00241 #ifdef USE_PARALLEL
00242    int comm_duped = 0;          /* Whether the MPI Communicator was duplicated */
00243    int info_duped = 0;          /* Whether the MPI Info object was duplicated */
00244 #else /* !USE_PARALLEL */
00245    int persist = 0; /* Should diskless try to persist its data into file?*/
00246 #endif
00247 
00248    assert(nc);
00249 
00250    if(cmode & NC_DISKLESS)
00251        flags = H5F_ACC_TRUNC;
00252    else if(cmode & NC_NOCLOBBER)
00253        flags = H5F_ACC_EXCL;
00254    else
00255        flags = H5F_ACC_TRUNC;
00256 
00257    LOG((3, "%s: path %s mode 0x%x", __func__, path, cmode));
00258    assert(nc && path);
00259 
00260    /* If this file already exists, and NC_NOCLOBBER is specified,
00261       return an error. */
00262    if (cmode & NC_DISKLESS) {
00263 #ifndef USE_PARALLEL
00264         if(cmode & NC_WRITE)
00265             persist = 1;
00266 #endif
00267    } else if ((cmode & NC_NOCLOBBER) && (fp = fopen(path, "r"))) {
00268       fclose(fp);
00269       return NC_EEXIST;
00270    }
00271    
00272    /* Add necessary structs to hold netcdf-4 file data. */
00273    if ((retval = nc4_nc4f_list_add(nc, path, (NC_WRITE | cmode))))
00274       BAIL(retval);
00275    nc4_info = NC4_DATA(nc);
00276    assert(nc4_info && nc4_info->root_grp);
00277 
00278 #if 0 /*def USE_PNETCDF*/
00279     if (cmode & NC_PNETCDF)
00280         return NC_NOERR;
00281 #endif
00282 
00283    /* Need this access plist to control how HDF5 handles open onjects
00284     * on file close. (Setting H5F_CLOSE_SEMI will cause H5Fclose to
00285     * fail if there are any open objects in the file. */
00286    if ((fapl_id = H5Pcreate(H5P_FILE_ACCESS)) < 0)
00287       BAIL(NC_EHDFERR);
00288 #ifdef EXTRA_TESTS
00289    num_plists++;
00290 #endif
00291 #ifdef EXTRA_TESTS
00292    if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_SEMI))
00293       BAIL(NC_EHDFERR);
00294 #else
00295    if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_STRONG))
00296       BAIL(NC_EHDFERR);
00297 #endif /* EXTRA_TESTS */
00298 
00299 #ifdef USE_PARALLEL
00300    /* If this is a parallel file create, set up the file creation
00301       property list. */
00302    if ((cmode & NC_MPIIO) || (cmode & NC_MPIPOSIX))
00303    {
00304       nc4_info->parallel = NC_TRUE;
00305       if (cmode & NC_MPIIO)  /* MPI/IO */
00306       {
00307          LOG((4, "creating parallel file with MPI/IO"));
00308          if (H5Pset_fapl_mpio(fapl_id, comm, info) < 0)
00309             BAIL(NC_EPARINIT);
00310       }
00311 #ifdef USE_PARALLEL_POSIX
00312       else /* MPI/POSIX */
00313       {
00314          LOG((4, "creating parallel file with MPI/posix"));
00315          if (H5Pset_fapl_mpiposix(fapl_id, comm, 0) < 0)
00316             BAIL(NC_EPARINIT);
00317       }
00318 #else /* USE_PARALLEL_POSIX */
00319       /* Should not happen! Code in NC4_create/NC4_open should alias the
00320        *        NC_MPIPOSIX flag to NC_MPIIO, if the MPI-POSIX VFD is not
00321        *        available in HDF5. -QAK
00322        */
00323       else /* MPI/POSIX */
00324          BAIL(NC_EPARINIT);
00325 #endif /* USE_PARALLEL_POSIX */
00326 
00327       /* Keep copies of the MPI Comm & Info objects */
00328       if (MPI_SUCCESS != MPI_Comm_dup(comm, &nc4_info->comm))
00329          BAIL(NC_EMPI);
00330       comm_duped++;
00331       if (MPI_INFO_NULL != info)
00332       {
00333          if (MPI_SUCCESS != MPI_Info_dup(info, &nc4_info->info))
00334             BAIL(NC_EMPI);
00335          info_duped++;
00336       }
00337       else
00338       {
00339          /* No dup, just copy it. */
00340          nc4_info->info = info;
00341       }
00342    }
00343 #else /* only set cache for non-parallel... */
00344    if(cmode & NC_DISKLESS) {
00345          if (H5Pset_fapl_core(fapl_id, 4096, persist))
00346             BAIL(NC_EDISKLESS);
00347    }
00348    if (H5Pset_cache(fapl_id, 0, nc4_chunk_cache_nelems, nc4_chunk_cache_size, 
00349                     nc4_chunk_cache_preemption) < 0)
00350       BAIL(NC_EHDFERR);
00351    LOG((4, "%s: set HDF raw chunk cache to size %d nelems %d preemption %f", 
00352         __func__, nc4_chunk_cache_size, nc4_chunk_cache_nelems, nc4_chunk_cache_preemption));
00353 #endif /* USE_PARALLEL */
00354    
00355    if (H5Pset_libver_bounds(fapl_id, H5F_LIBVER_LATEST, H5F_LIBVER_LATEST) < 0)
00356       BAIL(NC_EHDFERR);
00357 
00358    /* Create the property list. */
00359    if ((fcpl_id = H5Pcreate(H5P_FILE_CREATE)) < 0)
00360       BAIL(NC_EHDFERR);
00361 #ifdef EXTRA_TESTS
00362    num_plists++;
00363 #endif
00364 
00365    /* RJ: this suppose to be FALSE that is defined in H5 private.h as 0 */
00366    if (H5Pset_obj_track_times(fcpl_id,0)<0)
00367       BAIL(NC_EHDFERR);
00368 
00369    /* Set latest_format in access propertly list and
00370     * H5P_CRT_ORDER_TRACKED in the creation property list. This turns
00371     * on HDF5 creation ordering. */
00372    if (H5Pset_link_creation_order(fcpl_id, (H5P_CRT_ORDER_TRACKED |
00373                                             H5P_CRT_ORDER_INDEXED)) < 0)
00374       BAIL(NC_EHDFERR);
00375    if (H5Pset_attr_creation_order(fcpl_id, (H5P_CRT_ORDER_TRACKED |
00376                                             H5P_CRT_ORDER_INDEXED)) < 0)
00377       BAIL(NC_EHDFERR);
00378 
00379    /* Create the file. */
00380    if ((nc4_info->hdfid = H5Fcreate(path, flags, fcpl_id, fapl_id)) < 0) 
00381         /*Change the return error from NC_EFILEMETADATA to
00382           System error EACCES because that is the more likely problem */
00383       BAIL(EACCES);
00384 
00385    /* Open the root group. */
00386    if ((nc4_info->root_grp->hdf_grpid = H5Gopen2(nc4_info->hdfid, "/", 
00387                                                      H5P_DEFAULT)) < 0)
00388       BAIL(NC_EFILEMETA);
00389 
00390    /* Release the property lists. */
00391    if (H5Pclose(fapl_id) < 0 ||
00392        H5Pclose(fcpl_id) < 0)
00393       BAIL(NC_EHDFERR);
00394 #ifdef EXTRA_TESTS
00395    num_plists--;
00396    num_plists--;
00397 #endif
00398 
00399    /* Define mode gets turned on automatically on create. */
00400    nc4_info->flags |= NC_INDEF;
00401 
00402    return NC_NOERR;
00403 
00404 exit: /*failure exit*/
00405 #ifdef USE_PARALLEL
00406    if (comm_duped) MPI_Comm_free(&nc4_info->comm);
00407    if (info_duped) MPI_Info_free(&nc4_info->info);
00408 #endif
00409 #ifdef EXTRA_TESTS
00410    num_plists--;
00411 #endif
00412    if (fapl_id != H5P_DEFAULT) H5Pclose(fapl_id);
00413    if(!nc4_info) return retval;
00414    close_netcdf4_file(nc4_info,1); /* treat like abort */
00415    return retval;
00416 }
00417 
00433 int
00434 NC4_create(const char* path, int cmode, size_t initialsz, int basepe, 
00435            size_t *chunksizehintp, int use_parallel, void *mpidata,
00436            NC_Dispatch *dispatch, NC* nc_file)
00437 {
00438    MPI_Comm comm = MPI_COMM_WORLD; 
00439    MPI_Info info = MPI_INFO_NULL; 
00440    int res;
00441 
00442    assert(nc_file && path);
00443 
00444    LOG((1, "%s: path %s cmode 0x%x comm %d info %d",
00445         __func__, path, cmode, comm, info));
00446    
00447 #ifdef USE_PARALLEL
00448    if (mpidata) 
00449    { 
00450       comm = ((NC_MPI_INFO *)mpidata)->comm; 
00451       info = ((NC_MPI_INFO *)mpidata)->info;    
00452    }
00453 #endif /* USE_PARALLEL */
00454 
00455    /* If this is our first file, turn off HDF5 error messages. */
00456    if (virgin)
00457    {
00458       if (H5Eset_auto(NULL, NULL) < 0)
00459          LOG((0, "Couldn't turn off HDF5 error messages!"));
00460       LOG((1, "HDF5 error messages have been turned off."));
00461       virgin = 0;
00462    }
00463 
00464    /* Check the cmode for validity. */
00465    if (cmode & ~(NC_NOCLOBBER | NC_64BIT_OFFSET
00466                  | NC_NETCDF4 | NC_CLASSIC_MODEL
00467                  | NC_SHARE | NC_MPIIO | NC_MPIPOSIX | NC_LOCK | NC_PNETCDF
00468                  | NC_DISKLESS
00469                  | NC_WRITE /* to support diskless persistence */
00470                  )
00471        || (cmode & NC_MPIIO && cmode & NC_MPIPOSIX)
00472        || (cmode & NC_64BIT_OFFSET && cmode & NC_NETCDF4)
00473        || (cmode & (NC_MPIIO | NC_MPIPOSIX) && cmode & NC_DISKLESS)
00474       )
00475       return NC_EINVAL;
00476 
00477 #ifndef USE_PARALLEL_POSIX
00478 /* If the HDF5 library has been compiled without the MPI-POSIX VFD, alias
00479  *      the NC_MPIPOSIX flag to NC_MPIIO. -QAK
00480  */
00481    if(cmode & NC_MPIPOSIX)
00482    {
00483       cmode &= ~NC_MPIPOSIX;
00484       cmode |= NC_MPIIO;
00485    }
00486 #endif /* USE_PARALLEL_POSIX */
00487 
00488    cmode |= NC_NETCDF4;
00489 
00490    /* Apply default create format. */
00491    if (nc_get_default_format() == NC_FORMAT_64BIT)
00492       cmode |= NC_64BIT_OFFSET;
00493    else if (nc_get_default_format() == NC_FORMAT_NETCDF4_CLASSIC)
00494       cmode |= NC_CLASSIC_MODEL;
00495 
00496    LOG((2, "cmode after applying default format: 0x%x", cmode));
00497 
00498    nc_file->int_ncid = nc_file->ext_ncid;
00499    res = nc4_create_file(path, cmode, comm, info, nc_file);
00500 
00501 #if 0 /*def USE_PNETCDF*/
00502    if (cmode & NC_PNETCDF)
00503    {
00504       NC_HDF5_FILE_INFO_T* nc4_info;
00505       nc4_info = NC4_DATA(nc_file);
00506       assert(nc4_info);
00507 
00508       nc4_info->pnetcdf_file++;
00509       res = ncmpi_create(comm, path, cmode, info, &(nc_file->int_ncid));      
00510    }
00511 #endif /* USE_PNETCDF */
00512 
00513    return res;
00514 }
00515 
00516 /* This function is called by read_dataset when a dimension scale
00517  * dataset is encountered. It reads in the dimension data (creating a
00518  * new NC_DIM_INFO_T object), and also checks to see if this is a
00519  * dimension without a variable - that is, a coordinate dimension
00520  * which does not have any coordinate data. */
00521 static int
00522 read_scale(NC_GRP_INFO_T *grp, hid_t datasetid, const char *obj_name, 
00523         const H5G_stat_t *statbuf, hsize_t scale_size, hsize_t max_scale_size, 
00524         NC_DIM_INFO_T **dim)
00525 {
00526    NC_DIM_INFO_T *new_dim;              /* Dimension added to group */
00527    char dimscale_name_att[NC_MAX_NAME + 1] = "";    /* Dimscale name, for checking if dim without var */
00528    htri_t attr_exists = -1;             /* Flag indicating hidden attribute exists */
00529    hid_t attid = -1;                    /* ID of hidden attribute (to store dim ID) */
00530    int dimscale_created = 0;            /* Remember if a dimension was created (for error recovery) */
00531    int initial_grp_ndims = grp->ndims;  /* Retain for error recovery */
00532    short initial_next_dimid = grp->nc4_info->next_dimid;/* Retain for error recovery */
00533    int retval;
00534 
00535    /* Add a dimension for this scale. */
00536    if ((retval = nc4_dim_list_add(&grp->dim, &new_dim)))
00537       BAIL(retval);
00538    dimscale_created++;
00539 
00540    /* Does this dataset have a hidden attribute that tells us its
00541     * dimid? If so, read it. */
00542    if ((attr_exists = H5Aexists(datasetid, NC_DIMID_ATT_NAME)) < 0)
00543       BAIL(NC_EHDFERR);
00544    if (attr_exists)
00545    {
00546       if ((attid = H5Aopen_by_name(datasetid, ".", NC_DIMID_ATT_NAME, 
00547                                    H5P_DEFAULT, H5P_DEFAULT)) < 0)
00548          BAIL(NC_EHDFERR);
00549 
00550       if (H5Aread(attid, H5T_NATIVE_INT, &new_dim->dimid) < 0)
00551          BAIL(NC_EHDFERR);
00552 
00553       /* Check if scale's dimid should impact the group's next dimid */
00554       if (new_dim->dimid >= grp->nc4_info->next_dimid)
00555          grp->nc4_info->next_dimid = new_dim->dimid + 1;
00556    }
00557    else
00558    {
00559       /* Assign dimid */
00560       new_dim->dimid = grp->nc4_info->next_dimid++;
00561    }
00562 
00563    /* Increment number of dimensions. */
00564    grp->ndims++;
00565 
00566    if (!(new_dim->name = strdup(obj_name)))
00567       BAIL(NC_ENOMEM);
00568    if (SIZEOF_SIZE_T < 8 && scale_size > NC_MAX_UINT)
00569    {
00570       new_dim->len = NC_MAX_UINT;
00571       new_dim->too_long = NC_TRUE;
00572    }
00573    else
00574       new_dim->len = scale_size;
00575    new_dim->hdf5_objid.fileno[0] = statbuf->fileno[0];
00576    new_dim->hdf5_objid.fileno[1] = statbuf->fileno[1];
00577    new_dim->hdf5_objid.objno[0] = statbuf->objno[0];
00578    new_dim->hdf5_objid.objno[1] = statbuf->objno[1];
00579 
00580    /* If the dimscale has an unlimited dimension, then this dimension
00581     * is unlimited. */
00582    if (max_scale_size == H5S_UNLIMITED)
00583       new_dim->unlimited = NC_TRUE;
00584 
00585    /* If the scale name is set to DIM_WITHOUT_VARIABLE, then this is a
00586     * dimension, but not a variable. (If get_scale_name returns an
00587     * error, just move on, there's no NAME.) */
00588    if (H5DSget_scale_name(datasetid, dimscale_name_att, NC_MAX_NAME) >= 0)
00589    {
00590       if (!strncmp(dimscale_name_att, DIM_WITHOUT_VARIABLE, 
00591                    strlen(DIM_WITHOUT_VARIABLE)))
00592       {
00593          if (new_dim->unlimited)
00594          {
00595             size_t len = 0, *lenp = &len;
00596 
00597             if ((retval = nc4_find_dim_len(grp, new_dim->dimid, &lenp)))
00598                BAIL(retval);
00599             new_dim->len = *lenp;
00600          }
00601 
00602          /* Hold open the dataset, since the dimension doesn't have a coordinate variable */
00603          new_dim->hdf_dimscaleid = datasetid;
00604          H5Iinc_ref(new_dim->hdf_dimscaleid);        /* Increment number of objects using ID */
00605       }
00606    }
00607 
00608    /* Set the dimension created */
00609    *dim = new_dim;
00610 
00611 exit:
00612    /* Close the hidden attribute, if it was opened (error, or no error) */
00613    if (attid > 0 && H5Aclose(attid) < 0)
00614       BAIL2(NC_EHDFERR);
00615 
00616    /* On error, undo any dimscale creation */
00617    if (retval < 0 && dimscale_created)
00618    {
00619        /* Delete the dimension */
00620        if ((retval = nc4_dim_list_del(&grp->dim, new_dim)))
00621            BAIL2(retval);
00622 
00623       /* Reset the group's information */
00624       grp->ndims = initial_grp_ndims;
00625       grp->nc4_info->next_dimid = initial_next_dimid;
00626    }
00627 
00628    return retval;
00629 }
00630 
00631 /* This function reads the hacked in coordinates attribute I use for
00632  * multi-dimensional coordinates. */
00633 static int
00634 read_coord_dimids(NC_VAR_INFO_T *var)
00635 {
00636    hid_t coord_att_typeid = -1, coord_attid = -1, spaceid = -1;
00637    hssize_t npoints;
00638    int ret = 0;
00639 
00640    /* There is a hidden attribute telling us the ids of the
00641     * dimensions that apply to this multi-dimensional coordinate
00642     * variable. Read it. */
00643    if ((coord_attid = H5Aopen_name(var->hdf_datasetid, COORDINATES)) < 0) ret++;
00644    if (!ret && (coord_att_typeid = H5Aget_type(coord_attid)) < 0) ret++;
00645 
00646    /* How many dimensions are there? */
00647    if (!ret && (spaceid = H5Aget_space(coord_attid)) < 0) ret++;
00648 #ifdef EXTRA_TESTS
00649    num_spaces++;
00650 #endif
00651    if (!ret && (npoints = H5Sget_simple_extent_npoints(spaceid)) < 0) ret++;
00652 
00653    /* Check that the number of points is the same as the number of dimensions
00654     *   for the variable */
00655    if (!ret && npoints != var->ndims) ret++;
00656 
00657    if (!ret && H5Aread(coord_attid, coord_att_typeid, var->dimids) < 0) ret++;
00658    LOG((4, "dimscale %s is multidimensional and has coords", var->name));
00659    
00660    /* Set my HDF5 IDs free! */
00661    if (spaceid >= 0 && H5Sclose(spaceid) < 0) ret++;
00662 #ifdef EXTRA_TESTS
00663    num_spaces--;
00664 #endif
00665    if (coord_att_typeid >= 0 && H5Tclose(coord_att_typeid) < 0) ret++;
00666    if (coord_attid >= 0 && H5Aclose(coord_attid) < 0) ret++;
00667    return ret ? NC_EATTMETA : NC_NOERR;
00668 }
00669 
00670 /* This function is called when reading a file's metadata for each
00671  * dimension scale attached to a variable.*/
00672 static herr_t 
00673 dimscale_visitor(hid_t did, unsigned dim, hid_t dsid, 
00674                  void *dimscale_hdf5_objids)
00675 {
00676    H5G_stat_t statbuf;
00677 
00678    /* Get more info on the dimscale object.*/
00679    if (H5Gget_objinfo(dsid, ".", 1, &statbuf) < 0)
00680       return -1;
00681 
00682    /* Pass this information back to caller. */
00683    (*(HDF5_OBJID_T *)dimscale_hdf5_objids).fileno[0] = statbuf.fileno[0];
00684    (*(HDF5_OBJID_T *)dimscale_hdf5_objids).fileno[1] = statbuf.fileno[1];
00685    (*(HDF5_OBJID_T *)dimscale_hdf5_objids).objno[0] = statbuf.objno[0];
00686    (*(HDF5_OBJID_T *)dimscale_hdf5_objids).objno[1] = statbuf.objno[1];
00687    return 0;
00688 }
00689 
00690 /* Given an HDF5 type, set a pointer to netcdf type. */
00691 static int
00692 get_netcdf_type(NC_HDF5_FILE_INFO_T *h5, hid_t native_typeid, 
00693                 nc_type *xtype)
00694 {
00695    NC_TYPE_INFO_T *type;
00696    H5T_class_t class;
00697    htri_t is_str, equal = 0;
00698 
00699    assert(h5 && xtype);
00700 
00701    if ((class = H5Tget_class(native_typeid)) < 0)
00702       return NC_EHDFERR;
00703 
00704    /* H5Tequal doesn't work with H5T_C_S1 for some reason. But
00705     * H5Tget_class will return H5T_STRING if this is a string. */
00706    if (class == H5T_STRING)
00707    {
00708       if ((is_str = H5Tis_variable_str(native_typeid)) < 0)
00709          return NC_EHDFERR;
00710       if (is_str)
00711          *xtype = NC_STRING;
00712       else
00713          *xtype = NC_CHAR;
00714       return NC_NOERR;
00715    }
00716    else if (class == H5T_INTEGER || class == H5T_FLOAT)
00717    {
00718       /* For integers and floats, we don't have to worry about
00719        * endianness if we compare native types. */
00720       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_SCHAR)) < 0)
00721          return NC_EHDFERR;
00722       if (equal)
00723       {
00724          *xtype = NC_BYTE;
00725          return NC_NOERR;
00726       }
00727       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_SHORT)) < 0)
00728          return NC_EHDFERR;
00729       if (equal)
00730       {
00731          *xtype = NC_SHORT;
00732          return NC_NOERR;
00733       }
00734       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_INT)) < 0)
00735          return NC_EHDFERR;
00736       if (equal)
00737       {
00738          *xtype = NC_INT;
00739          return NC_NOERR;
00740       }
00741       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_FLOAT)) < 0)
00742          return NC_EHDFERR;
00743       if (equal)
00744       {
00745          *xtype = NC_FLOAT;
00746          return NC_NOERR;
00747       }
00748       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_DOUBLE)) < 0)
00749          return NC_EHDFERR;
00750       if (equal)
00751       {
00752          *xtype = NC_DOUBLE;
00753          return NC_NOERR;
00754       }
00755       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_UCHAR)) < 0)
00756          return NC_EHDFERR;
00757       if (equal)
00758       {
00759          *xtype = NC_UBYTE;
00760          return NC_NOERR;
00761       }
00762       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_USHORT)) < 0)
00763          return NC_EHDFERR;
00764       if (equal)
00765       {
00766          *xtype = NC_USHORT;
00767          return NC_NOERR;
00768       }
00769       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_UINT)) < 0)
00770          return NC_EHDFERR;
00771       if (equal)
00772       {
00773          *xtype = NC_UINT;
00774          return NC_NOERR;
00775       }
00776       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_LLONG)) < 0)
00777          return NC_EHDFERR;
00778       if (equal)
00779       {
00780          *xtype = NC_INT64;
00781          return NC_NOERR;
00782       }
00783       if ((equal = H5Tequal(native_typeid, H5T_NATIVE_ULLONG)) < 0)
00784          return NC_EHDFERR;
00785       if (equal)
00786       {
00787          *xtype = NC_UINT64;
00788          return NC_NOERR;
00789       }
00790    }
00791 
00792    /* Maybe we already know about this type. */
00793    if (!equal)
00794       if((type = nc4_rec_find_hdf_type(h5->root_grp, native_typeid)))
00795       {
00796          *xtype = type->nc_typeid;
00797          return NC_NOERR;
00798       }
00799    
00800    *xtype = NC_NAT;
00801    return NC_EBADTYPID;
00802 }
00803 
00804 /* Given an HDF5 type, set a pointer to netcdf type_info struct,
00805  * either an existing one (for user-defined types) or a newly created
00806  * one. */
00807 static int
00808 get_type_info2(NC_HDF5_FILE_INFO_T *h5, hid_t datasetid,
00809                NC_TYPE_INFO_T **type_info)
00810 {
00811    htri_t is_str, equal = 0;
00812    H5T_class_t class;
00813    hid_t native_typeid, hdf_typeid;
00814    H5T_order_t order;
00815    int t;
00816 
00817    assert(h5 && type_info);
00818 
00819    /* Because these N5T_NATIVE_* constants are actually function calls
00820     * (!) in H5Tpublic.h, I can't initialize this array in the usual
00821     * way, because at least some C compilers (like Irix) complain
00822     * about calling functions when defining constants. So I have to do
00823     * it like this. Note that there's no native types for char or
00824     * string. Those are handled later. */
00825    if (!h5_native_type_constant_g[1])
00826    {
00827       h5_native_type_constant_g[1] = H5T_NATIVE_SCHAR;
00828       h5_native_type_constant_g[2] = H5T_NATIVE_SHORT;
00829       h5_native_type_constant_g[3] = H5T_NATIVE_INT;
00830       h5_native_type_constant_g[4] = H5T_NATIVE_FLOAT;
00831       h5_native_type_constant_g[5] = H5T_NATIVE_DOUBLE;
00832       h5_native_type_constant_g[6] = H5T_NATIVE_UCHAR;
00833       h5_native_type_constant_g[7] = H5T_NATIVE_USHORT;
00834       h5_native_type_constant_g[8] = H5T_NATIVE_UINT;
00835       h5_native_type_constant_g[9] = H5T_NATIVE_LLONG;
00836       h5_native_type_constant_g[10] = H5T_NATIVE_ULLONG;
00837    }
00838    
00839    /* Get the HDF5 typeid - we'll need it later. */
00840    if ((hdf_typeid = H5Dget_type(datasetid)) < 0)
00841       return NC_EHDFERR;
00842 
00843    /* Get the native typeid. Will be equivalent to hdf_typeid when
00844     * creating but not necessarily when reading, a variable. */
00845    if ((native_typeid = H5Tget_native_type(hdf_typeid, H5T_DIR_DEFAULT)) < 0) 
00846       return NC_EHDFERR;
00847 
00848    /* Is this type an integer, string, compound, or what? */
00849    if ((class = H5Tget_class(native_typeid)) < 0)
00850       return NC_EHDFERR;
00851 
00852    /* Is this an atomic type? */
00853    if (class == H5T_STRING || class == H5T_INTEGER || class == H5T_FLOAT)
00854    {
00855       /* Allocate a phony NC_TYPE_INFO_T struct to hold type info. */
00856       if (!(*type_info = calloc(1, sizeof(NC_TYPE_INFO_T))))
00857          return NC_ENOMEM;
00858 
00859       /* H5Tequal doesn't work with H5T_C_S1 for some reason. But
00860        * H5Tget_class will return H5T_STRING if this is a string. */
00861       if (class == H5T_STRING)
00862       {
00863          if ((is_str = H5Tis_variable_str(native_typeid)) < 0)
00864             return NC_EHDFERR;
00865          /* Make sure fixed-len strings will work like variable-len strings */
00866          if (is_str || H5Tget_size(hdf_typeid) > 1)
00867          {
00868             /* Set a class for the type */
00869             t = NUM_TYPES - 1;
00870             (*type_info)->nc_type_class = NC_STRING;
00871          }
00872          else
00873          {
00874             /* Set a class for the type */
00875             t = 0;
00876             (*type_info)->nc_type_class = NC_CHAR;
00877          }
00878       }
00879       else if (class == H5T_INTEGER || class == H5T_FLOAT)
00880       {
00881          for (t = 1; t < NUM_TYPES - 1; t++)
00882          {
00883             if ((equal = H5Tequal(native_typeid, h5_native_type_constant_g[t])) < 0)
00884                return NC_EHDFERR;
00885             if (equal)
00886                break;
00887          }
00888 
00889          /* Find out about endianness. */
00890          if (class == H5T_INTEGER)
00891          {
00892             if ((order = H5Tget_order(hdf_typeid)) < 0) 
00893                return NC_EHDFERR;
00894 
00895             /* Copy this into the type_info struct. */
00896             if (order == H5T_ORDER_LE)
00897                (*type_info)->endianness = NC_ENDIAN_LITTLE;
00898             else if (order == H5T_ORDER_BE)
00899                (*type_info)->endianness = NC_ENDIAN_BIG;
00900             else /* don't support H5T_ORDER_VAX, H5T_ORDER_MIXED, H5T_ORDER_NONE */
00901                 return NC_EBADTYPE;
00902 
00903             /* Set a class for the type */
00904             /* (Note use of 'NC_INT' for all integer types) */
00905             (*type_info)->nc_type_class = NC_INT;
00906          }
00907          else
00908          {
00909             /* Set a class for the type */
00910             /* (Note use of 'NC_FLOAT' for all floating-point types) */
00911             (*type_info)->nc_type_class = NC_FLOAT;
00912          }
00913       }
00914       (*type_info)->nc_typeid = nc_type_constant_g[t];
00915       (*type_info)->size = nc_type_size_g[t];
00916       if (!((*type_info)->name = strdup(nc_type_name_g[t])))
00917          return NC_ENOMEM;
00918       (*type_info)->hdf_typeid = hdf_typeid;
00919       (*type_info)->native_hdf_typeid = native_typeid;
00920       return NC_NOERR;
00921    }
00922    else
00923    {
00924       NC_TYPE_INFO_T *type;
00925 
00926       /* This is a user-defined type. */
00927       if((type = nc4_rec_find_hdf_type(h5->root_grp, native_typeid)))
00928          *type_info = type;
00929 
00930       /* The type entry in the array of user-defined types already has
00931        * an open data typeid (and native typeid), so close the ones we
00932        * opened above. */
00933       if (H5Tclose(native_typeid) < 0) 
00934          return NC_EHDFERR;
00935       if (H5Tclose(hdf_typeid) < 0) 
00936          return NC_EHDFERR;
00937 
00938       if (type)
00939          return NC_NOERR;
00940    }
00941 
00942    return NC_EBADTYPID;
00943 }
00944 
00945 /* Read an attribute. */
00946 static int 
00947 read_hdf5_att(NC_GRP_INFO_T *grp, hid_t attid, NC_ATT_INFO_T *att)
00948 {
00949    hid_t spaceid = 0, file_typeid = 0;
00950    hsize_t dims[1] = {0}; /* netcdf attributes always 1-D. */
00951    int retval = NC_NOERR;
00952    size_t type_size;
00953    int att_ndims;
00954    hssize_t att_npoints;
00955    H5T_class_t att_class;      
00956    int fixed_len_string = 0;
00957    size_t fixed_size = 0;
00958 
00959    assert(att->name);
00960    LOG((5, "%s: att->attnum %d att->name %s att->nc_typeid %d att->len %d",
00961         __func__, att->attnum, att->name, (int)att->nc_typeid, att->len));
00962 
00963    /* Get type of attribute in file. */
00964    if ((file_typeid = H5Aget_type(attid)) < 0)
00965       return NC_EATTMETA;
00966    if ((att->native_hdf_typeid = H5Tget_native_type(file_typeid, H5T_DIR_DEFAULT)) < 0) 
00967       BAIL(NC_EHDFERR);
00968    if ((att_class = H5Tget_class(att->native_hdf_typeid)) < 0)
00969       BAIL(NC_EATTMETA);
00970    if (att_class == H5T_STRING && !H5Tis_variable_str(att->native_hdf_typeid))
00971    {
00972       fixed_len_string++;
00973       if (!(fixed_size = H5Tget_size(att->native_hdf_typeid)))
00974          BAIL(NC_EATTMETA);
00975    }
00976    if ((retval = get_netcdf_type(grp->nc4_info, att->native_hdf_typeid, &(att->nc_typeid))))
00977       BAIL(retval);
00978 
00979 
00980    /* Get len. */
00981    if ((spaceid = H5Aget_space(attid)) < 0)
00982       BAIL(NC_EATTMETA); 
00983 #ifdef EXTRA_TESTS
00984    num_spaces++;
00985 #endif
00986    if ((att_ndims = H5Sget_simple_extent_ndims(spaceid)) < 0)
00987       BAIL(NC_EATTMETA);
00988    if ((att_npoints = H5Sget_simple_extent_npoints(spaceid)) < 0)
00989       BAIL(NC_EATTMETA);
00990 
00991    /* If both att_ndims and att_npoints are zero, then this is a
00992     * zero length att. */
00993    if (att_ndims == 0 && att_npoints == 0)
00994       dims[0] = 0;
00995    else if (att->nc_typeid == NC_STRING)
00996        dims[0] = att_npoints;
00997    else if (att->nc_typeid == NC_CHAR)
00998    {
00999       /* NC_CHAR attributes are written as a scalar in HDF5, of type
01000        * H5T_C_S1, of variable length. */
01001       if (att_ndims == 0) 
01002       {
01003          if (!(dims[0] = H5Tget_size(file_typeid)))
01004             BAIL(NC_EATTMETA);
01005       }
01006       else
01007       {
01008          /* This is really a string type! */
01009          att->nc_typeid = NC_STRING;
01010          dims[0] = att_npoints;
01011       }
01012    } 
01013    else
01014    {
01015       H5S_class_t space_class;
01016 
01017       /* All netcdf attributes are scalar or 1-D only. */
01018       if (att_ndims > 1)
01019          BAIL(NC_EATTMETA);
01020 
01021       /* Check class of HDF5 dataspace */
01022       if ((space_class = H5Sget_simple_extent_type(spaceid)) < 0)
01023          BAIL(NC_EATTMETA);
01024 
01025       /* Check for NULL HDF5 dataspace class (should be weeded out earlier) */
01026       if (H5S_NULL == space_class)
01027          BAIL(NC_EATTMETA);
01028 
01029       /* check for SCALAR HDF5 dataspace class */
01030       if (H5S_SCALAR == space_class)
01031           dims[0] = 1;
01032       else /* Must be "simple" dataspace */
01033       {
01034           /* Read the size of this attribute. */
01035           if (H5Sget_simple_extent_dims(spaceid, dims, NULL) < 0)
01036              BAIL(NC_EATTMETA);
01037       }
01038    }
01039       
01040    /* Tell the user what the length if this attribute is. */
01041    att->len = dims[0];
01042 
01043    /* Allocate some memory if the len is not zero, and read the
01044       attribute. */
01045    if (dims[0])
01046    {
01047       if ((retval = nc4_get_typelen_mem(grp->nc4_info, att->nc_typeid, 0,
01048                                         &type_size)))
01049          return retval;
01050       if (att_class == H5T_VLEN)
01051       {
01052          if (!(att->vldata = malloc((unsigned int)(att->len * sizeof(hvl_t)))))
01053             BAIL(NC_ENOMEM);
01054          if (H5Aread(attid, att->native_hdf_typeid, att->vldata) < 0)
01055             BAIL(NC_EATTMETA);
01056       }
01057       else if (att->nc_typeid == NC_STRING)
01058       {
01059          if (!(att->stdata = calloc(att->len, sizeof(char *))))
01060             BAIL(NC_ENOMEM);
01061          /* For a fixed length HDF5 string, the read requires
01062           * contiguous memory. Meanwhile, the netCDF API requires that
01063           * nc_free_string be called on string arrays, which would not
01064           * work if one contiguous memory block were used. So here I
01065           * convert the contiguous block of strings into an array of
01066           * malloced strings (each string with its own malloc). Then I
01067           * copy the data and free the contiguous memory. This
01068           * involves copying the data, which is bad, but this only
01069           * occurs for fixed length string attributes, and presumably
01070           * these are small. (And netCDF-4 does not create them - it
01071           * always uses variable length strings. */
01072          if (fixed_len_string)
01073          {
01074             int i;
01075             char *contig_buf, *cur;
01076 
01077             /* Alloc space for the contiguous memory read. */
01078             if (!(contig_buf = malloc(att->len * fixed_size * sizeof(char))))
01079                BAIL(NC_ENOMEM);
01080 
01081             /* Read the fixed-len strings as one big block. */
01082             if (H5Aread(attid, att->native_hdf_typeid, contig_buf) < 0)
01083                BAIL(NC_EATTMETA);
01084             
01085             /* Copy strings, one at a time, into their new home. Alloc
01086                space for each string. The user will later free this
01087                space with nc_free_string. */
01088             cur = contig_buf;
01089             for (i = 0; i < att->len; i++)
01090             {
01091                if (!(att->stdata[i] = malloc(fixed_size)))
01092                   BAIL(NC_ENOMEM);
01093                strncpy(att->stdata[i], cur, fixed_size);
01094                cur += fixed_size;
01095             }
01096             
01097             /* Free contiguous memory buffer. */
01098             free(contig_buf);
01099          }
01100          else
01101          {
01102             /* Read variable-length string atts. */
01103             if (H5Aread(attid, att->native_hdf_typeid, att->stdata) < 0)
01104                BAIL(NC_EATTMETA);
01105          }
01106       }
01107       else
01108       {
01109          if (!(att->data = malloc((unsigned int)(att->len * type_size))))
01110             BAIL(NC_ENOMEM);
01111          if (H5Aread(attid, att->native_hdf_typeid, att->data) < 0)
01112             BAIL(NC_EATTMETA);
01113       }
01114    }
01115 
01116    if (H5Tclose(file_typeid) < 0)
01117       BAIL(NC_EHDFERR);
01118    if (H5Sclose(spaceid) < 0)
01119       return NC_EHDFERR;
01120 #ifdef EXTRA_TESTS
01121    num_spaces--;
01122 #endif
01123    
01124    return NC_NOERR;
01125 
01126   exit:
01127    if (H5Tclose(file_typeid) < 0)
01128       BAIL2(NC_EHDFERR);
01129    if (spaceid > 0 && H5Sclose(spaceid) < 0)
01130       BAIL2(NC_EHDFERR);
01131 #ifdef EXTRA_TESTS
01132    num_spaces--;
01133 #endif
01134    return retval;
01135 }
01136 
01137 /* Read information about a user defined type from the HDF5 file, and
01138  * stash it in the group's list of types. */
01139 static int
01140 read_type(NC_GRP_INFO_T *grp, hid_t hdf_typeid, char *type_name)
01141 {
01142    NC_TYPE_INFO_T *type;
01143    H5T_class_t class;
01144    hid_t native_typeid;
01145    size_t type_size;
01146    int retval = NC_NOERR;
01147 
01148    assert(grp && type_name);
01149 
01150    LOG((4, "%s: type_name %s grp->name %s", __func__, type_name, grp->name));
01151 
01152    /* What is the native type for this platform? */
01153    if ((native_typeid = H5Tget_native_type(hdf_typeid, H5T_DIR_DEFAULT)) < 0) 
01154       return NC_EHDFERR;
01155    
01156    /* What is the size of this type on this platform. */
01157    if (!(type_size = H5Tget_size(native_typeid)))
01158       return NC_EHDFERR;
01159    LOG((5, "type_size %d", type_size));
01160 
01161    /* Add to the list for this new type, and get a local pointer to it. */
01162    if ((retval = nc4_type_list_add(grp, type_size, type_name, &type)))
01163       return retval;
01164 
01165    /* Remember common info about this type. */
01166    type->committed = NC_TRUE;
01167    type->hdf_typeid = hdf_typeid;
01168    H5Iinc_ref(type->hdf_typeid);        /* Increment number of objects using ID */
01169    type->native_hdf_typeid = native_typeid;
01170 
01171    /* What is the class of this type, compound, vlen, etc. */
01172    if ((class = H5Tget_class(hdf_typeid)) < 0)
01173       return NC_EHDFERR;
01174    switch (class)
01175    {
01176       case H5T_STRING:
01177          type->nc_type_class = NC_STRING;
01178          break;
01179 
01180       case H5T_COMPOUND:
01181          {
01182             int nmembers;
01183             unsigned int m;
01184             char* member_name = NULL;
01185 #ifdef JNA
01186             char jna[1001];
01187 #endif  
01188 
01189             type->nc_type_class = NC_COMPOUND;
01190 
01191             if ((nmembers = H5Tget_nmembers(hdf_typeid)) < 0)
01192                return NC_EHDFERR;
01193             LOG((5, "compound type has %d members", nmembers));
01194             for (m = 0; m < nmembers; m++)
01195             {
01196                hid_t member_hdf_typeid;
01197                hid_t member_native_typeid;
01198                size_t member_offset;
01199                H5T_class_t mem_class;
01200                nc_type member_xtype;
01201 
01202                retval = NC_NOERR;
01203 
01204                /* Get the typeid and native typeid of this member of the
01205                 * compound type. */
01206                if ((member_hdf_typeid = H5Tget_member_type(type->native_hdf_typeid, m)) < 0)
01207                   return NC_EHDFERR;
01208 
01209                if ((member_native_typeid = H5Tget_native_type(member_hdf_typeid, H5T_DIR_DEFAULT)) < 0) 
01210                   return NC_EHDFERR;
01211 
01212                /* Get the name of the member.*/
01213                member_name = H5Tget_member_name(type->native_hdf_typeid, m);
01214                if (!member_name || strlen(member_name) > NC_MAX_NAME) {
01215                   retval = NC_EBADNAME;
01216                   break;
01217                }
01218 #ifdef JNA
01219                else {
01220                 strncpy(jna,member_name,1000);
01221                 member_name = jna;      
01222                }
01223 #endif
01224 
01225                /* Offset in bytes on *this* platform. */
01226                member_offset = H5Tget_member_offset(type->native_hdf_typeid, m);
01227 
01228                /* Get dimensional data if this member is an array of something. */
01229                if ((mem_class = H5Tget_class(member_hdf_typeid)) < 0)
01230                   return NC_EHDFERR;
01231                if (mem_class == H5T_ARRAY)
01232                {
01233                   int ndims, dim_size[NC_MAX_VAR_DIMS];
01234                   hsize_t dims[NC_MAX_VAR_DIMS];
01235                   int d;
01236 
01237                   if ((ndims = H5Tget_array_ndims(member_hdf_typeid)) < 0) {
01238                      retval = NC_EHDFERR;
01239                      break;
01240                   }
01241                   if (H5Tget_array_dims(member_hdf_typeid, dims, NULL) != ndims) {
01242                      retval = NC_EHDFERR;
01243                      break;
01244                   }
01245                   for (d = 0; d < ndims; d++)
01246                      dim_size[d] = dims[d];
01247 
01248                   /* What is the netCDF typeid of this member? */
01249                   if ((retval = get_netcdf_type(grp->nc4_info, H5Tget_super(member_hdf_typeid), 
01250                                                 &member_xtype)))
01251                      break;
01252 
01253                   /* Add this member to our list of fields in this compound type. */
01254                   if ((retval = nc4_field_list_add(&type->u.c.field, type->u.c.num_fields++, member_name, 
01255                                                    member_offset, H5Tget_super(member_hdf_typeid), 
01256                                                    H5Tget_super(member_native_typeid), 
01257                                                    member_xtype, ndims, dim_size)))
01258                      break;
01259                }
01260                else
01261                {
01262                   /* What is the netCDF typeid of this member? */
01263                   if ((retval = get_netcdf_type(grp->nc4_info, member_native_typeid, 
01264                                                 &member_xtype)))
01265                      break;
01266 
01267                   /* Add this member to our list of fields in this compound type. */
01268                   if ((retval = nc4_field_list_add(&type->u.c.field, type->u.c.num_fields++, member_name, 
01269                                                    member_offset, member_hdf_typeid, member_native_typeid, 
01270                                                    member_xtype, 0, NULL)))
01271                      break;
01272                } 
01273                
01274 #ifndef JNA
01275                /* Free the member name (which HDF5 allocated for us). */
01276                if(member_name != NULL) free(member_name); 
01277 #endif         
01278                member_name = NULL;
01279             }
01280 #ifndef JNA
01281             if(member_name != NULL)
01282                 free(member_name);
01283 #endif
01284             if(retval) /* error exit from loop */
01285                 return retval;
01286          }
01287          break;
01288 
01289       case H5T_VLEN:
01290          {
01291             htri_t ret;
01292 
01293             /* For conveninence we allow user to pass vlens of strings
01294              * with null terminated strings. This means strings are
01295              * treated slightly differently by the API, although they are
01296              * really just VLENs of characters. */
01297             if ((ret = H5Tis_variable_str(hdf_typeid)) < 0)
01298                return NC_EHDFERR;
01299             if (ret)
01300                type->nc_type_class = NC_STRING;
01301             else
01302             {
01303                hid_t base_hdf_typeid;
01304                nc_type base_nc_type = NC_NAT;
01305 
01306                type->nc_type_class = NC_VLEN;
01307 
01308                /* Find the base type of this vlen (i.e. what is this a
01309                 * vlen of?) */
01310                if (!(base_hdf_typeid = H5Tget_super(native_typeid)))
01311                   return NC_EHDFERR;
01312 
01313                /* What size is this type? */
01314                if (!(type_size = H5Tget_size(base_hdf_typeid)))
01315                   return NC_EHDFERR;
01316 
01317                /* What is the netcdf corresponding type. */
01318                if ((retval = get_netcdf_type(grp->nc4_info, base_hdf_typeid, 
01319                                              &base_nc_type)))
01320                   return retval;
01321                LOG((5, "base_hdf_typeid 0x%x type_size %d base_nc_type %d", 
01322                     base_hdf_typeid, type_size, base_nc_type));
01323 
01324                /* Remember the base types for this vlen */
01325                type->u.v.base_nc_typeid = base_nc_type;
01326                type->u.v.base_hdf_typeid = base_hdf_typeid;
01327             }
01328          }
01329          break;
01330 
01331       case H5T_OPAQUE:
01332          type->nc_type_class = NC_OPAQUE;
01333          break;
01334 
01335       case H5T_ENUM:
01336          {
01337             hid_t base_hdf_typeid;
01338             nc_type base_nc_type = NC_NAT;
01339             void *value;
01340             int i;
01341             char *member_name = NULL;
01342 #ifdef JNA
01343             char jna[1001];
01344 #endif  
01345 
01346             type->nc_type_class = NC_ENUM;
01347 
01348             /* Find the base type of this enum (i.e. what is this a
01349              * enum of?) */
01350             if (!(base_hdf_typeid = H5Tget_super(hdf_typeid)))
01351                return NC_EHDFERR;
01352             /* What size is this type? */
01353             if (!(type_size = H5Tget_size(base_hdf_typeid)))
01354                return NC_EHDFERR;
01355             /* What is the netcdf corresponding type. */
01356             if ((retval = get_netcdf_type(grp->nc4_info, base_hdf_typeid, 
01357                                           &base_nc_type)))
01358                return retval;
01359             LOG((5, "base_hdf_typeid 0x%x type_size %d base_nc_type %d", 
01360                  base_hdf_typeid, type_size, base_nc_type));
01361 
01362             /* Remember the base types for this enum */
01363             type->u.e.base_nc_typeid = base_nc_type;
01364             type->u.e.base_hdf_typeid = base_hdf_typeid;
01365 
01366             /* Find out how many member are in the enum. */
01367             if ((type->u.e.num_members = H5Tget_nmembers(hdf_typeid)) < 0) 
01368                return NC_EHDFERR;
01369 
01370             /* Allocate space for one value. */
01371             if (!(value = calloc(1, type_size)))
01372                return NC_ENOMEM;
01373 
01374             /* Read each name and value defined in the enum. */
01375             for (i = 0; i < type->u.e.num_members; i++)
01376             {
01377                retval = NC_NOERR;
01378                /* Get the name and value from HDF5. */
01379                if (!(member_name = H5Tget_member_name(hdf_typeid, i)))
01380                {
01381                   retval = NC_EHDFERR;
01382                   break;                  
01383                }
01384 #ifdef JNA
01385                 strncpy(jna,member_name,1000);
01386                 member_name = jna;      
01387 #endif
01388 
01389                if (strlen(member_name) > NC_MAX_NAME)
01390                {
01391                   retval = NC_EBADNAME;
01392                   break;
01393                }
01394                if (H5Tget_member_value(hdf_typeid, i, value) < 0) 
01395                {
01396                   retval = NC_EHDFERR;
01397                   break;
01398                }
01399 
01400                /* Insert new field into this type's list of fields. */
01401                if ((retval = nc4_enum_member_add(&type->u.e.enum_member, type->size, 
01402                                                  member_name, value)))
01403                {
01404                   break;
01405                }
01406 
01407 #ifndef JNA
01408                /* Free the member name (which HDF5 allocated for us). */
01409                if(member_name != NULL) free(member_name); 
01410 #endif         
01411                member_name = NULL;
01412             }
01413 
01414 #ifndef JNA
01415             if(member_name != NULL)
01416                 free(member_name);
01417 #endif
01418             if(value) free(value);
01419             if(retval) /* error exit from loop */
01420                 return retval;
01421          }
01422          break;
01423 
01424       default:
01425          LOG((0, "unknown class"));
01426          return NC_EBADCLASS;
01427    }
01428    return retval;
01429 }
01430 
01431 /* This function is called by read_dataset, (which is called by
01432  * nc4_rec_read_metadata) when a netCDF variable is found in the
01433  * file. This function reads in all the metadata about the var,
01434  * including the attributes. */
01435 static int
01436 read_var(NC_GRP_INFO_T *grp, hid_t datasetid, const char *obj_name, 
01437          size_t ndims, NC_DIM_INFO_T *dim)
01438 {
01439    NC_VAR_INFO_T *var = NULL;
01440    hid_t access_pid = 0;
01441    int incr_id_rc = 0;          /* Whether the dataset ID's ref count has been incremented */
01442    int natts, a, d;
01443 
01444    NC_ATT_INFO_T *att;
01445    hid_t attid = 0;
01446    char att_name[NC_MAX_HDF5_NAME + 1];
01447 
01448 #define CD_NELEMS_ZLIB 1
01449 #define CD_NELEMS_SZIP 4
01450    H5Z_filter_t filter;
01451    int num_filters;
01452    unsigned int cd_values[CD_NELEMS_SZIP];
01453    size_t cd_nelems = CD_NELEMS_SZIP;
01454    hid_t propid = 0;
01455    H5D_fill_value_t fill_status;
01456    H5D_layout_t layout;
01457    hsize_t chunksize[NC_MAX_VAR_DIMS] = {0};
01458    int retval = NC_NOERR;
01459    double rdcc_w0;
01460    int f;
01461 
01462    assert(obj_name && grp);
01463    LOG((4, "%s: obj_name %s", __func__, obj_name));
01464 
01465    /* Add a variable to the end of the group's var list. */
01466    if ((retval = nc4_var_list_add(&grp->var, &var)))
01467       BAIL(retval);
01468    
01469    /* Fill in what we already know. */
01470    var->hdf_datasetid = datasetid;
01471    H5Iinc_ref(var->hdf_datasetid);      /* Increment number of objects using ID */
01472    incr_id_rc++;                        /* Indicate that we've incremented the ref. count (for errors) */
01473    var->varid = grp->nvars++;
01474    var->created = NC_TRUE;
01475    var->ndims = ndims;
01476 
01477    /* We need some room to store information about dimensions for this
01478     * var. */
01479    if (var->ndims)
01480    {
01481       if (!(var->dim = calloc(var->ndims, sizeof(NC_DIM_INFO_T *))))
01482          BAIL(NC_ENOMEM);
01483       if (!(var->dimids = calloc(var->ndims, sizeof(int))))
01484          BAIL(NC_ENOMEM);
01485    }
01486 
01487    /* Get the current chunk cache settings. */
01488    if ((access_pid = H5Dget_access_plist(datasetid)) < 0)
01489       BAIL(NC_EVARMETA);
01490 #ifdef EXTRA_TESTS
01491    num_plists++;
01492 #endif
01493 
01494    /* Learn about current chunk cache settings. */
01495    if ((H5Pget_chunk_cache(access_pid, &(var->chunk_cache_nelems), 
01496                            &(var->chunk_cache_size), &rdcc_w0)) < 0)
01497       BAIL(NC_EHDFERR);
01498    var->chunk_cache_preemption = rdcc_w0;
01499 
01500    /* Check for a weird case: a non-coordinate variable that has the
01501     * same name as a dimension. It's legal in netcdf, and requires
01502     * that the HDF5 dataset name be changed. */
01503    if (strlen(obj_name) > strlen(NON_COORD_PREPEND) &&
01504          !strncmp(obj_name, NON_COORD_PREPEND, strlen(NON_COORD_PREPEND)))
01505    {
01506       /* Allocate space for the name. */
01507       if (!(var->name = malloc(((strlen(obj_name) - strlen(NON_COORD_PREPEND))+ 1) * sizeof(char))))
01508          BAIL(NC_ENOMEM);
01509 
01510       strcpy(var->name, &obj_name[strlen(NON_COORD_PREPEND)]);
01511 
01512       /* Allocate space for the HDF5 name. */
01513       if (!(var->hdf5_name = malloc((strlen(obj_name) + 1) * sizeof(char))))
01514          BAIL(NC_ENOMEM);
01515 
01516       strcpy(var->hdf5_name, obj_name);
01517    }
01518    else
01519    {
01520       /* Allocate space for the name. */
01521       if (!(var->name = malloc((strlen(obj_name) + 1) * sizeof(char))))
01522          BAIL(NC_ENOMEM);
01523 
01524       strcpy(var->name, obj_name);
01525    }
01526 
01527    /* Find out what filters are applied to this HDF5 dataset,
01528     * fletcher32, deflate, and/or shuffle. All other filters are
01529     * ignored. */
01530    if ((propid = H5Dget_create_plist(datasetid)) < 0) 
01531       BAIL(NC_EHDFERR);
01532 #ifdef EXTRA_TESTS
01533    num_plists++;
01534 #endif /* EXTRA_TESTS */
01535 
01536    /* Get the chunking info for non-scalar vars. */
01537    if ((layout = H5Pget_layout(propid)) < -1)
01538       BAIL(NC_EHDFERR);
01539    if (layout == H5D_CHUNKED)
01540    {
01541       if (H5Pget_chunk(propid, NC_MAX_VAR_DIMS, chunksize) < 0)
01542          BAIL(NC_EHDFERR);
01543       if (!(var->chunksizes = malloc(var->ndims * sizeof(size_t))))
01544          BAIL(NC_ENOMEM);
01545       for (d = 0; d < var->ndims; d++)
01546          var->chunksizes[d] = chunksize[d];
01547    }
01548    else if (layout == H5D_CONTIGUOUS || layout == H5D_COMPACT)
01549       var->contiguous = NC_TRUE;
01550 
01551    /* The possible values of filter (which is just an int) can be
01552     * found in H5Zpublic.h. */
01553    if ((num_filters = H5Pget_nfilters(propid)) < 0) 
01554       BAIL(NC_EHDFERR);
01555    for (f = 0; f < num_filters; f++)
01556    {
01557       if ((filter = H5Pget_filter2(propid, f, NULL, &cd_nelems, 
01558                                    cd_values, 0, NULL, NULL)) < 0)
01559          BAIL(NC_EHDFERR);
01560       switch (filter)
01561       {
01562          case H5Z_FILTER_SHUFFLE:
01563             var->shuffle = NC_TRUE;
01564             break;
01565 
01566          case H5Z_FILTER_FLETCHER32:
01567             var->fletcher32 = NC_TRUE;
01568             break;
01569 
01570          case H5Z_FILTER_DEFLATE:
01571             var->deflate = NC_TRUE;
01572             if (cd_nelems != CD_NELEMS_ZLIB || cd_values[0] > MAX_DEFLATE_LEVEL)
01573                BAIL(NC_EHDFERR);
01574             var->deflate_level = cd_values[0];
01575             break;
01576 
01577          case H5Z_FILTER_SZIP:
01578             var->szip = NC_TRUE;
01579             if (cd_nelems != CD_NELEMS_SZIP)
01580                BAIL(NC_EHDFERR);
01581             var->options_mask = cd_values[0];
01582             var->pixels_per_block = cd_values[1];
01583             break;
01584 
01585          default:
01586             LOG((1, "Yikes! Unknown filter type found on dataset!"));
01587             break;
01588       }
01589    }
01590                
01591    /* Learn all about the type of this variable. */
01592    if ((retval = get_type_info2(grp->nc4_info, datasetid, 
01593                                 &var->type_info)))
01594       BAIL(retval);
01595 
01596    /* Indicate that the variable has a pointer to the type */
01597    var->type_info->rc++;
01598 
01599    /* Is there a fill value associated with this dataset? */
01600    if (H5Pfill_value_defined(propid, &fill_status) < 0)
01601       BAIL(NC_EHDFERR);
01602 
01603    /* Get the fill value, if there is one defined. */
01604    if (fill_status == H5D_FILL_VALUE_USER_DEFINED)
01605    {
01606       /* Allocate space to hold the fill value. */
01607       if (!var->fill_value)
01608       {
01609          if (var->type_info->nc_type_class == NC_VLEN)
01610          {
01611             if (!(var->fill_value = malloc(sizeof(nc_vlen_t))))
01612                BAIL(NC_ENOMEM);
01613          }
01614          else if (var->type_info->nc_type_class == NC_STRING)
01615          {
01616             if (!(var->fill_value = malloc(sizeof(char *))))
01617                BAIL(NC_ENOMEM);
01618          }
01619          else
01620          {
01621             assert(var->type_info->size);
01622             if (!(var->fill_value = malloc(var->type_info->size)))
01623                BAIL(NC_ENOMEM);
01624          }
01625       }
01626       
01627       /* Get the fill value from the HDF5 property lust. */
01628       if (H5Pget_fill_value(propid, var->type_info->native_hdf_typeid, 
01629                             var->fill_value) < 0)
01630          BAIL(NC_EHDFERR);
01631    }
01632    else
01633       var->no_fill = NC_TRUE;
01634 
01635    /* If it's a scale, mark it as such. */
01636    if (dim)
01637    {
01638       assert(ndims);
01639       var->dimscale = NC_TRUE;
01640       if (var->ndims > 1)
01641       {
01642          if ((retval = read_coord_dimids(var)))
01643             BAIL(retval);
01644       }
01645       else
01646       {
01647          /* sanity check */
01648          assert(0 == strcmp(var->name, dim->name));
01649 
01650          var->dimids[0] = dim->dimid;
01651          var->dim[0] = dim;
01652       }
01653       dim->coord_var = var;
01654    }
01655    /* If this is not a scale, but has scales, iterate
01656     * through them. (i.e. this is a variable that is not a
01657     * coordinate variable) */
01658    else
01659    {
01660       int num_scales = 0;            
01661 
01662       /* Find out how many scales are attached to this
01663        * dataset. H5DSget_num_scales returns an error if there are no
01664        * scales, so convert a negative return value to zero. */
01665       num_scales = H5DSget_num_scales(datasetid, 0);
01666       if (num_scales < 0)
01667          num_scales = 0;
01668 
01669       if (num_scales && ndims)
01670       {
01671          /* Allocate space to remember whether the dimscale has been attached
01672           * for each dimension. */
01673          if (NULL == (var->dimscale_attached = calloc(ndims, sizeof(nc_bool_t))))
01674             BAIL(NC_ENOMEM);       
01675 
01676          /* Store id information allowing us to match hdf5
01677           * dimscales to netcdf dimensions. */
01678          if (NULL == (var->dimscale_hdf5_objids = malloc(ndims * sizeof(struct hdf5_objid))))
01679             BAIL(NC_ENOMEM);
01680          for (d = 0; d < var->ndims; d++)
01681          {
01682             if (H5DSiterate_scales(var->hdf_datasetid, d, NULL, dimscale_visitor,
01683                                    &(var->dimscale_hdf5_objids[d])) < 0)
01684                BAIL(NC_EHDFERR);
01685             var->dimscale_attached[d] = NC_TRUE;
01686          }
01687       }
01688    }
01689         
01690    /* Now read all the attributes of this variable, ignoring the
01691       ones that hold HDF5 dimension scale information. */
01692    if ((natts = H5Aget_num_attrs(var->hdf_datasetid)) < 0)
01693       BAIL(NC_EATTMETA);
01694    for (a = 0; a < natts; a++)
01695    {
01696       /* Close the attribute and try to move on with our
01697        * lives. Like bits through the network port, so
01698        * flows the Days of Our Lives! */
01699       if (attid && H5Aclose(attid) < 0)
01700          BAIL(NC_EHDFERR);
01701 
01702       /* Open the att and get its name. */
01703       if ((attid = H5Aopen_idx(var->hdf_datasetid, (unsigned int)a)) < 0)
01704          BAIL(NC_EATTMETA);
01705       if (H5Aget_name(attid, NC_MAX_HDF5_NAME, att_name) < 0)
01706          BAIL(NC_EATTMETA);
01707       LOG((4, "%s:: a %d att_name %s", __func__, a, att_name));
01708 
01709       /* Should we ignore this attribute? */    
01710       if (strcmp(att_name, REFERENCE_LIST) &&
01711           strcmp(att_name, CLASS) &&
01712           strcmp(att_name, DIMENSION_LIST) &&
01713           strcmp(att_name, NAME) &&
01714           strcmp(att_name, COORDINATES) &&
01715           strcmp(att_name, NC_DIMID_ATT_NAME))
01716       {
01717          /* Add to the end of the list of atts for this var. */
01718          if ((retval = nc4_att_list_add(&var->att, &att)))
01719             BAIL(retval);
01720          
01721          /* Fill in the information we know. */
01722          att->attnum = var->natts++;
01723          if (!(att->name = strdup(att_name)))
01724             BAIL(NC_ENOMEM);
01725          
01726          /* Read the rest of the info about the att,
01727           * including its values. */
01728          if ((retval = read_hdf5_att(grp, attid, att)))
01729          {
01730             if (NC_EBADTYPID == retval)
01731             {
01732                 if ((retval = nc4_att_list_del(&var->att, att)))
01733                     BAIL(retval);
01734                 continue;
01735             }
01736             else
01737                 BAIL(retval);
01738          }
01739          
01740          att->created = NC_TRUE;
01741       } /* endif not HDF5 att */
01742    } /* next attribute */
01743 
01744    /* Is this a deflated variable with a chunksize greater than the
01745     * current cache size? */
01746    if ((retval = nc4_adjust_var_cache(grp, var)))
01747       BAIL(retval);
01748 
01749 exit:
01750    if (retval)
01751    {
01752        if (incr_id_rc && H5Idec_ref(datasetid) < 0)
01753           BAIL2(NC_EHDFERR);
01754        if (var && nc4_var_list_del(&grp->var, var))
01755           BAIL2(NC_EHDFERR);
01756    }
01757    if (access_pid && H5Pclose(access_pid) < 0)
01758       BAIL2(NC_EHDFERR);
01759 #ifdef EXTRA_TESTS
01760    num_plists--;
01761 #endif
01762    if (propid > 0 && H5Pclose(propid) < 0)
01763       BAIL2(NC_EHDFERR);
01764 #ifdef EXTRA_TESTS
01765    num_plists--;
01766 #endif
01767    if (attid > 0 && H5Aclose(attid) < 0)
01768       BAIL2(NC_EHDFERR);
01769    return retval;
01770 }
01771 
01772 /* This function is called by nc4_rec_read_metadata to read all the
01773  * group level attributes (the NC_GLOBAL atts for this group). */
01774 static int
01775 read_grp_atts(NC_GRP_INFO_T *grp)
01776 {
01777    hid_t attid = 0;
01778    hsize_t num_obj, i;
01779    NC_ATT_INFO_T *att;
01780    NC_TYPE_INFO_T *type;
01781    char obj_name[NC_MAX_HDF5_NAME + 1];
01782    int max_len;
01783    int retval = NC_NOERR;
01784 
01785    num_obj = H5Aget_num_attrs(grp->hdf_grpid);
01786    for (i = 0; i < num_obj; i++)
01787    {
01788       /* Close an attribute from previous loop iteration */
01789       /* (Should be from 'continue' statement, below) */
01790       if (attid && H5Aclose(attid) < 0)
01791          BAIL(NC_EHDFERR);
01792 
01793       if ((attid = H5Aopen_idx(grp->hdf_grpid, (unsigned int)i)) < 0)
01794          BAIL(NC_EATTMETA);
01795       if (H5Aget_name(attid, NC_MAX_NAME + 1, obj_name) < 0)
01796          BAIL(NC_EATTMETA);
01797       LOG((3, "reading attribute of _netCDF group, named %s", obj_name));
01798 
01799       /* This may be an attribute telling us that strict netcdf-3
01800        * rules are in effect. If so, we will make note of the fact,
01801        * but not add this attribute to the metadata. It's not a user
01802        * attribute, but an internal netcdf-4 one. */
01803       if (!strcmp(obj_name, NC3_STRICT_ATT_NAME))
01804          grp->nc4_info->cmode |= NC_CLASSIC_MODEL;
01805       else
01806       {
01807          /* Add an att struct at the end of the list, and then go to it. */
01808          if ((retval = nc4_att_list_add(&grp->att, &att)))
01809             BAIL(retval);
01810 
01811          /* Add the info about this attribute. */
01812          max_len = strlen(obj_name) > NC_MAX_NAME ? NC_MAX_NAME : strlen(obj_name);
01813          if (!(att->name = malloc((max_len + 1) * sizeof(char))))
01814             BAIL(NC_ENOMEM);
01815          strncpy(att->name, obj_name, max_len);
01816          att->name[max_len] = 0;
01817          att->attnum = grp->natts++;
01818          if ((retval = read_hdf5_att(grp, attid, att)))
01819          {
01820             if (NC_EBADTYPID == retval)
01821             {
01822                if ((retval = nc4_att_list_del(&grp->att, att)))
01823                   BAIL(retval);
01824                continue;
01825             }
01826             else
01827                BAIL(retval);
01828          }
01829          att->created = NC_TRUE;
01830          if ((retval = nc4_find_type(grp->nc4_info, att->nc_typeid, &type)))
01831             BAIL(retval);
01832       }
01833    }
01834 
01835   exit:
01836    if (attid > 0 && H5Aclose(attid) < 0)
01837       BAIL2(NC_EHDFERR);
01838    return retval;
01839 }
01840 
01841 /* This function is called when nc4_rec_read_metadata encounters an HDF5
01842  * dataset when reading a file. */
01843 static int
01844 read_dataset(NC_GRP_INFO_T *grp, hid_t datasetid, const char *obj_name,
01845     const H5G_stat_t *statbuf)
01846 {
01847    NC_DIM_INFO_T *dim = NULL;   /* Dimension created for scales */
01848    hid_t spaceid = 0;
01849    int ndims;
01850    int is_scale = 0;
01851    int retval = NC_NOERR;
01852 
01853    /* Get the dimension information for this dataset. */
01854    if ((spaceid = H5Dget_space(datasetid)) < 0)
01855       BAIL(NC_EHDFERR);
01856 #ifdef EXTRA_TESTS
01857    num_spaces++;
01858 #endif
01859    if ((ndims = H5Sget_simple_extent_ndims(spaceid)) < 0)
01860       BAIL(NC_EHDFERR);
01861 
01862    /* Is this a dimscale? */
01863    if ((is_scale = H5DSis_scale(datasetid)) < 0)
01864       BAIL(NC_EHDFERR);
01865    if (is_scale)
01866    {
01867       hsize_t dims[H5S_MAX_RANK];
01868       hsize_t max_dims[H5S_MAX_RANK];
01869 
01870       /* Query the scale's size & max. size */
01871       if (H5Sget_simple_extent_dims(spaceid, dims, max_dims) < 0)
01872          BAIL(NC_EHDFERR);
01873 
01874       /* Read the scale information. */
01875       if ((retval = read_scale(grp, datasetid, obj_name, statbuf, dims[0],
01876                                max_dims[0], &dim)))
01877          BAIL(retval);
01878    }
01879 
01880    /* Add a var to the linked list, and get its metadata,
01881     * unless this is one of those funny dimscales that are a
01882     * dimension in netCDF but not a variable. (Spooky!) */
01883    if (NULL == dim || (dim && !dim->hdf_dimscaleid))
01884       if ((retval = read_var(grp, datasetid, obj_name, ndims, dim)))
01885          BAIL(retval);
01886    
01887 exit: 
01888    if (spaceid && H5Sclose(spaceid) <0)
01889       BAIL2(retval);
01890 #ifdef EXTRA_TESTS
01891    num_spaces--;
01892 #endif
01893 
01894    return retval;
01895 }
01896 
01897 static int
01898 nc4_rec_read_metadata_cb_list_add(NC4_rec_read_metadata_obj_info_t **head,
01899                   NC4_rec_read_metadata_obj_info_t **tail,
01900                   const NC4_rec_read_metadata_obj_info_t *oinfo)
01901 {
01902    NC4_rec_read_metadata_obj_info_t *new_oinfo;    /* Pointer to info for object */
01903 
01904    /* Allocate memory for the object's info */
01905    if (!(new_oinfo = calloc(1, sizeof(*new_oinfo))))
01906       return NC_ENOMEM;
01907 
01908    /* Make a copy of the object's info */
01909    memcpy(new_oinfo, oinfo, sizeof(*oinfo));
01910 
01911    if (*tail)
01912    {
01913        assert(*head);
01914        (*tail)->next = new_oinfo;
01915        *tail = new_oinfo;
01916    }
01917    else
01918    {
01919        assert(NULL == *head);
01920        *head = *tail = new_oinfo;
01921    }
01922 
01923    return (NC_NOERR);
01924 }
01925 
01926 static int
01927 nc4_rec_read_metadata_cb(hid_t grpid, const char *name, const H5L_info_t *info,
01928                       void *_op_data)
01929 {
01930     NC4_rec_read_metadata_ud_t *udata = (NC4_rec_read_metadata_ud_t *)_op_data; /* Pointer to user data for callback */
01931     NC4_rec_read_metadata_obj_info_t oinfo;    /* Pointer to info for object */
01932     int retval = H5_ITER_CONT;
01933 
01934    /* Reset the memory for the object's info */
01935    memset(&oinfo, 0, sizeof(oinfo));
01936 
01937    /* Open this critter. */
01938    if ((oinfo.oid = H5Oopen(grpid, name, H5P_DEFAULT)) < 0) 
01939       BAIL(H5_ITER_ERROR);
01940           
01941    /* Get info about the object.*/
01942    if (H5Gget_objinfo(oinfo.oid, ".", 1, &oinfo.statbuf) < 0)
01943       BAIL(H5_ITER_ERROR);
01944 
01945    strncpy(oinfo.oname, name, NC_MAX_NAME);
01946            
01947    /* Add object to list, for later */
01948    switch(oinfo.statbuf.type)
01949    {
01950       case H5G_GROUP:
01951          LOG((3, "found group %s", oinfo.oname));
01952 
01953          /* Defer descending into child group immediately, so that the types
01954           *     in the current group can be processed and be ready for use by
01955           *     vars in the child group(s).
01956           */
01957          if (nc4_rec_read_metadata_cb_list_add(&udata->grps_head, &udata->grps_tail, &oinfo))
01958              BAIL(H5_ITER_ERROR);
01959          break;
01960 
01961       case H5G_DATASET:
01962          LOG((3, "found dataset %s", oinfo.oname));
01963 
01964          /* Learn all about this dataset, which may be a dimscale
01965           * (i.e. dimension metadata), or real data. */
01966          if ((retval = read_dataset(udata->grp, oinfo.oid, oinfo.oname, &oinfo.statbuf)))
01967          {
01968             /* Allow NC_EBADTYPID to transparently skip over datasets
01969              *  which have a datatype that netCDF-4 doesn't undertand
01970              *  (currently), but break out of iteration for other
01971              *  errors.
01972              */
01973             if(NC_EBADTYPID != retval)
01974                BAIL(H5_ITER_ERROR);
01975             else
01976                retval = H5_ITER_CONT;
01977          }
01978 
01979          /* Close the object */
01980          if (H5Oclose(oinfo.oid) < 0)
01981             BAIL(H5_ITER_ERROR);
01982          break;
01983 
01984       case H5G_TYPE:
01985          LOG((3, "found datatype %s", oinfo.oname));
01986 
01987          /* Process the named datatype */
01988          if (read_type(udata->grp, oinfo.oid, oinfo.oname))
01989             BAIL(H5_ITER_ERROR);
01990 
01991          /* Close the object */
01992          if (H5Oclose(oinfo.oid) < 0)
01993             BAIL(H5_ITER_ERROR);
01994          break;
01995 
01996       default:
01997          LOG((0, "Unknown object class %d in %s!", oinfo.statbuf.type, __func__));
01998          BAIL(H5_ITER_ERROR);
01999     }
02000 
02001 exit:
02002    if (retval)
02003    {
02004       if (oinfo.oid > 0 && H5Oclose(oinfo.oid) < 0)
02005          BAIL2(H5_ITER_ERROR);
02006    }
02007 
02008    return (retval);
02009 }
02010 
02011 /* This is the main function to recursively read all the metadata for the file. */
02012 /* The links in the 'grp' are iterated over and added to the file's metadata
02013  *      information.  Note that child groups are not immediately processed, but
02014  *      are deferred until all the other links in the group are handled (so that
02015  *      vars in the child groups are guaranteed to have types that they use in
02016  *      a parent group in place).
02017  */
02018 static int
02019 nc4_rec_read_metadata(NC_GRP_INFO_T *grp)
02020 {
02021     NC4_rec_read_metadata_ud_t udata;   /* User data for iteration */
02022     NC4_rec_read_metadata_obj_info_t *oinfo;    /* Pointer to info for object */
02023     hsize_t idx=0;
02024     hid_t pid = 0;
02025     unsigned crt_order_flags = 0;
02026     H5_index_t iter_index;
02027     int retval = NC_NOERR; /* everything worked! */
02028 
02029     assert(grp && grp->name);
02030     LOG((3, "%s: grp->name %s", __func__, grp->name));
02031 
02032     /* Portably initialize user data for later */
02033     memset(&udata, 0, sizeof(udata));
02034 
02035     /* Open this HDF5 group and retain its grpid. It will remain open
02036      * with HDF5 until this file is nc_closed. */
02037     if (!grp->hdf_grpid)
02038     {
02039         if (grp->parent)
02040         {
02041             if ((grp->hdf_grpid = H5Gopen2(grp->parent->hdf_grpid, 
02042                                         grp->name, H5P_DEFAULT)) < 0)
02043                 BAIL(NC_EHDFERR);
02044         }
02045         else
02046         {
02047             if ((grp->hdf_grpid = H5Gopen2(grp->nc4_info->hdfid, 
02048                                            "/", H5P_DEFAULT)) < 0)
02049                 BAIL(NC_EHDFERR);
02050         }
02051     }
02052     assert(grp->hdf_grpid > 0);
02053 
02054     /* Get the group creation flags, to check for creation ordering */
02055     pid = H5Gget_create_plist(grp->hdf_grpid);
02056     H5Pget_link_creation_order(pid, &crt_order_flags); 
02057     if (H5Pclose(pid) < 0)
02058         BAIL(NC_EHDFERR);
02059         
02060     /* Set the iteration index to use */
02061     if (crt_order_flags & H5P_CRT_ORDER_TRACKED)
02062         iter_index = H5_INDEX_CRT_ORDER;
02063     else
02064     {
02065         NC_HDF5_FILE_INFO_T *h5 = grp->nc4_info;
02066 
02067         /* Without creation ordering, file must be read-only. */
02068         if (!h5->no_write)
02069             BAIL(NC_ECANTWRITE);
02070 
02071         iter_index = H5_INDEX_NAME;
02072     }
02073 
02074     /* Set user data for iteration */
02075     udata.grp = grp;
02076 
02077     /* Iterate over links in this group, building lists for the types,
02078      *  datasets and groups encountered
02079      */
02080     if (H5Literate(grp->hdf_grpid, iter_index, H5_ITER_INC, &idx,
02081             nc4_rec_read_metadata_cb, (void *)&udata) < 0)
02082         BAIL(NC_EHDFERR);
02083 
02084     /* Process the child groups found */
02085     /* (Deferred until now, so that the types in the current group get
02086      *  processed and are available for vars in the child group(s).)
02087      */
02088     for (oinfo = udata.grps_head; oinfo; oinfo = udata.grps_head)
02089     {
02090         NC_GRP_INFO_T *child_grp;
02091         NC_HDF5_FILE_INFO_T *h5 = grp->nc4_info;
02092 
02093         /* Add group to file's hierarchy */
02094         if ((retval = nc4_grp_list_add(&(grp->children), h5->next_nc_grpid++, 
02095                         grp, grp->nc4_info->controller, oinfo->oname, &child_grp)))
02096             BAIL(retval);
02097 
02098         /* Recursively read the child group's metadata */
02099         if ((retval = nc4_rec_read_metadata(child_grp)))
02100             BAIL(retval);
02101 
02102         /* Close the object */
02103         if (H5Oclose(oinfo->oid) < 0)
02104             BAIL(NC_EHDFERR);
02105 
02106         /* Advance to next node, free current node */
02107         udata.grps_head = oinfo->next;
02108         free(oinfo);
02109     }
02110 
02111     /* Scan the group for global (i.e. group-level) attributes. */
02112     if ((retval = read_grp_atts(grp)))
02113         BAIL(retval);
02114     
02115 exit:
02116     /* Clean up local information on error, if anything remains */
02117     if (retval)
02118     {
02119         for (oinfo = udata.grps_head; oinfo; oinfo = udata.grps_head)
02120         {
02121             /* Close the object */
02122             if (H5Oclose(oinfo->oid) < 0)
02123                 BAIL2(NC_EHDFERR);
02124 
02125             /* Advance to next node, free current node */
02126             udata.grps_head = oinfo->next;
02127             free(oinfo);
02128         }
02129     }
02130 
02131     return retval;
02132 }
02133 
02134 /* Open a netcdf-4 file. Things have already been kicked off in
02135  * ncfunc.c in nc_open, but here the netCDF-4 part of opening a file
02136  * is handled. */
02137 static int
02138 nc4_open_file(const char *path, int mode, MPI_Comm comm,
02139               MPI_Info info, NC *nc)
02140 {
02141    hid_t fapl_id = H5P_DEFAULT;
02142    unsigned flags = (mode & NC_WRITE) ? 
02143       H5F_ACC_RDWR : H5F_ACC_RDONLY;
02144    int retval;
02145    NC_HDF5_FILE_INFO_T* nc4_info = NULL;
02146 #ifdef USE_PARALLEL
02147    int comm_duped = 0;          /* Whether the MPI Communicator was duplicated */
02148    int info_duped = 0;          /* Whether the MPI Info object was duplicated */
02149 #endif /* !USE_PARALLEL */
02150 
02151    LOG((3, "%s: path %s mode %d", __func__, path, mode));
02152    assert(path && nc);
02153    /* Stop diskless open in its tracks */
02154    if(mode & NC_DISKLESS)
02155         return NC_EDISKLESS;
02156 
02157    /* Add necessary structs to hold netcdf-4 file data. */
02158    if ((retval = nc4_nc4f_list_add(nc, path, mode)))
02159       BAIL(retval);
02160    nc4_info = NC4_DATA(nc);
02161    assert(nc4_info && nc4_info->root_grp);
02162    
02163    /* Need this access plist to control how HDF5 handles open onjects
02164     * on file close. (Setting H5F_CLOSE_SEMI will cause H5Fclose to
02165     * fail if there are any open objects in the file. */
02166    if ((fapl_id = H5Pcreate(H5P_FILE_ACCESS)) < 0)
02167       BAIL(NC_EHDFERR);
02168 #ifdef EXTRA_TESTS
02169    num_plists++;
02170 #endif      
02171 #ifdef EXTRA_TESTS
02172    if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_SEMI)) 
02173       BAIL(NC_EHDFERR);
02174 #else
02175    if (H5Pset_fclose_degree(fapl_id, H5F_CLOSE_STRONG))
02176       BAIL(NC_EHDFERR);
02177 #endif
02178 
02179 #ifdef USE_PARALLEL
02180    /* If this is a parallel file create, set up the file creation
02181       property list. */
02182    if (mode & NC_MPIIO || mode & NC_MPIPOSIX)
02183    {
02184       nc4_info->parallel = NC_TRUE;
02185       if (mode & NC_MPIIO)  /* MPI/IO */
02186       {
02187          LOG((4, "opening parallel file with MPI/IO"));
02188          if (H5Pset_fapl_mpio(fapl_id, comm, info) < 0)
02189             BAIL(NC_EPARINIT);
02190       }
02191 #ifdef USE_PARALLEL_POSIX
02192       else /* MPI/POSIX */
02193       {
02194          LOG((4, "opening parallel file with MPI/posix"));
02195          if (H5Pset_fapl_mpiposix(fapl_id, comm, 0) < 0)
02196             BAIL(NC_EPARINIT);
02197       }
02198 #else /* USE_PARALLEL_POSIX */
02199       /* Should not happen! Code in NC4_create/NC4_open should alias the
02200        *        NC_MPIPOSIX flag to NC_MPIIO, if the MPI-POSIX VFD is not
02201        *        available in HDF5. -QAK
02202        */
02203       else /* MPI/POSIX */
02204          BAIL(NC_EPARINIT);
02205 #endif /* USE_PARALLEL_POSIX */
02206 
02207       /* Keep copies of the MPI Comm & Info objects */
02208       if (MPI_SUCCESS != MPI_Comm_dup(comm, &nc4_info->comm))
02209          BAIL(NC_EMPI);
02210       comm_duped++;
02211       if (MPI_INFO_NULL != info)
02212       {
02213          if (MPI_SUCCESS != MPI_Info_dup(info, &nc4_info->info))
02214             BAIL(NC_EMPI);
02215          info_duped++;
02216       }
02217       else
02218       {
02219          /* No dup, just copy it. */
02220          nc4_info->info = info;
02221       }
02222    }
02223 #else /* only set cache for non-parallel. */
02224    if (H5Pset_cache(fapl_id, 0, nc4_chunk_cache_nelems, nc4_chunk_cache_size, 
02225                     nc4_chunk_cache_preemption) < 0)
02226       BAIL(NC_EHDFERR);
02227    LOG((4, "%s: set HDF raw chunk cache to size %d nelems %d preemption %f", 
02228         __func__, nc4_chunk_cache_size, nc4_chunk_cache_nelems, nc4_chunk_cache_preemption));
02229 #endif /* USE_PARALLEL */
02230    
02231    /* The NetCDF-3.x prototype contains an mode option NC_SHARE for
02232       multiple processes accessing the dataset concurrently.  As there
02233       is no HDF5 equivalent, NC_SHARE is treated as NC_NOWRITE. */
02234    if ((nc4_info->hdfid = H5Fopen(path, flags, fapl_id)) < 0)
02235       BAIL(NC_EHDFERR);
02236 
02237    /* Does the mode specify that this file is read-only? */
02238    if ((mode & NC_WRITE) == 0)
02239       nc4_info->no_write = NC_TRUE;
02240 
02241    /* Now read in all the metadata. Some types and dimscale
02242     * information may be difficult to resolve here, if, for example, a
02243     * dataset of user-defined type is encountered before the
02244     * definition of that type. */
02245    if ((retval = nc4_rec_read_metadata(nc4_info->root_grp)))
02246       BAIL(retval);
02247 
02248    /* Now figure out which netCDF dims are indicated by the dimscale
02249     * information. */
02250    if ((retval = nc4_rec_match_dimscales(nc4_info->root_grp)))
02251       BAIL(retval);
02252 
02253 #ifdef LOGGING
02254    /* This will print out the names, types, lens, etc of the vars and
02255       atts in the file, if the logging level is 2 or greater. */ 
02256    log_metadata_nc(nc);
02257 #endif
02258 
02259    /* Close the property list. */ 
02260    if (H5Pclose(fapl_id) < 0)
02261       BAIL(NC_EHDFERR);
02262 #ifdef EXTRA_TESTS
02263    num_plists--;
02264 #endif
02265 
02266    return NC_NOERR;
02267 
02268 exit:
02269 #ifdef USE_PARALLEL
02270    if (comm_duped) MPI_Comm_free(&nc4_info->comm);
02271    if (info_duped) MPI_Info_free(&nc4_info->info);
02272 #endif
02273 #ifdef EXTRA_TESTS
02274    num_plists--;
02275 #endif
02276    if (fapl_id != H5P_DEFAULT) H5Pclose(fapl_id);
02277    if (!nc4_info) return retval;
02278    close_netcdf4_file(nc4_info,1); /*  treat like abort*/
02279    return retval;
02280 }
02281 
02282 /* Given an HDF4 type, set a pointer to netcdf type. */
02283 #ifdef USE_HDF4   
02284 static int
02285 get_netcdf_type_from_hdf4(NC_HDF5_FILE_INFO_T *h5, int32 hdf4_typeid, 
02286                           nc_type *xtype, NC_TYPE_INFO_T *type_info)
02287 {
02288    int t;
02289    assert(h5 && xtype);
02290 
02291    switch(hdf4_typeid)
02292    {
02293       case DFNT_CHAR:
02294          *xtype = NC_CHAR;
02295          t = 0;
02296          break;
02297       case DFNT_UCHAR:
02298       case DFNT_UINT8:
02299          *xtype = NC_UBYTE;
02300          t = 6;
02301          break;
02302       case DFNT_INT8:
02303          *xtype = NC_BYTE;
02304          t = 1;
02305          break;
02306       case DFNT_INT16:
02307          *xtype = NC_SHORT;
02308          t = 2;
02309          break;
02310       case DFNT_UINT16:
02311          *xtype = NC_USHORT;
02312          t = 7;
02313          break;
02314       case DFNT_INT32:
02315          *xtype = NC_INT;
02316          t = 3;
02317          break;
02318       case DFNT_UINT32:
02319          *xtype = NC_UINT;
02320          t = 8;
02321          break;
02322       case DFNT_FLOAT32:
02323          *xtype = NC_FLOAT;
02324          t = 4;
02325          break;
02326       case DFNT_FLOAT64:
02327          *xtype = NC_DOUBLE;
02328          t = 5;
02329          break;
02330       default:
02331          *xtype = NC_NAT;
02332          return NC_EBADTYPID;
02333    }
02334 
02335    if (type_info)
02336    {
02337       if (hdf4_typeid == DFNT_FLOAT32)
02338          type_info->nc_type_class = NC_FLOAT;
02339       else if (hdf4_typeid == DFNT_FLOAT64)
02340          type_info->nc_type_class = NC_DOUBLE;
02341       else if (hdf4_typeid == DFNT_CHAR)
02342          type_info->nc_type_class = NC_STRING;
02343       else
02344          type_info->nc_type_class = NC_INT;
02345       type_info->endianness = NC_ENDIAN_BIG;
02346       type_info->nc_typeid = *xtype;
02347       type_info->size = nc_type_size_g[t];
02348       if (!(type_info->name = strdup(nc_type_name_g[t])))
02349          return NC_ENOMEM;
02350    }
02351 
02352    return NC_NOERR;
02353 }
02354 
02355 /* Open a HDF4 file. Things have already been kicked off in nc_open,
02356  * but here the netCDF-4 part of opening a file is handled. */
02357 static int
02358 nc4_open_hdf4_file(const char *path, int mode, NC *nc)
02359 {
02360    NC_HDF5_FILE_INFO_T *h5;
02361    NC_GRP_INFO_T *grp;
02362    NC_ATT_INFO_T *att;
02363    int32 num_datasets, num_gatts;
02364    int32 rank;
02365    int v, d, a;
02366    int retval;
02367    NC_HDF5_FILE_INFO_T* nc4_info = NULL;
02368 
02369    LOG((3, "%s: path %s mode %d", __func__, path, mode));
02370    assert(path && nc);
02371 
02372    /* Must be read-only access to hdf4 files. */
02373    if (mode & NC_WRITE)
02374       return NC_EINVAL;
02375 
02376    /* Add necessary structs to hold netcdf-4 file data. */
02377    if ((retval = nc4_nc4f_list_add(nc, path, mode)))
02378       return retval;
02379    nc4_info = NC4_DATA(nc);
02380    assert(nc4_info && nc4_info->root_grp);
02381    h5 = nc4_info;
02382    h5->hdf4 = NC_TRUE;
02383    grp = h5->root_grp;
02384    h5->no_write = NC_TRUE;
02385 
02386    /* Open the file and initialize SD interface. */
02387    if ((h5->sdid = SDstart(path, DFACC_READ)) == FAIL)
02388       return NC_EHDFERR;
02389 
02390    /* Learn how many datasets and global atts we have. */
02391    if (SDfileinfo(h5->sdid, &num_datasets, &num_gatts))
02392       return NC_EHDFERR;
02393 
02394    /* Read the atts. */
02395    for (a = 0; a < num_gatts; a++)
02396    {
02397       int32 att_data_type, att_count;
02398       size_t att_type_size;
02399 
02400       /* Add to the end of the list of atts for this var. */
02401       if ((retval = nc4_att_list_add(&h5->root_grp->att, &att)))
02402          return retval;
02403       att->attnum = grp->natts++;
02404       att->created = NC_TRUE;
02405 
02406       /* Learn about this attribute. */
02407       if (!(att->name = malloc(NC_MAX_HDF4_NAME * sizeof(char))))
02408          return NC_ENOMEM;
02409       if (SDattrinfo(h5->sdid, a, att->name, &att_data_type, &att_count)) 
02410          return NC_EATTMETA;
02411       if ((retval = get_netcdf_type_from_hdf4(h5, att_data_type, 
02412                                               &att->nc_typeid, NULL)))
02413          return retval;
02414       att->len = att_count;
02415 
02416       /* Allocate memory to hold the data. */
02417       if ((retval = nc4_get_typelen_mem(h5, att->nc_typeid, 0, &att_type_size)))
02418          return retval;
02419       if (!(att->data = malloc(att_type_size * att->len)))
02420          return NC_ENOMEM;
02421 
02422       /* Read the data. */
02423       if (SDreadattr(h5->sdid, a, att->data)) 
02424          return NC_EHDFERR;
02425    }
02426 
02427    /* Read each dataset. */
02428    for (v = 0; v < num_datasets; v++)
02429    {
02430       NC_VAR_INFO_T *var;
02431       int32 data_type, num_atts;
02432       /* Problem: Number of dims is returned by the call that requires
02433          a pre-allocated array, 'dimsize'. 
02434        From SDS_SD website: 
02435        http://www.hdfgroup.org/training/HDFtraining/UsersGuide/SDS_SD.fm3.html 
02436        The maximum rank is 32, or MAX_VAR_DIMS (as defined in netcdf.h).
02437        
02438        int32 dimsize[MAX_VAR_DIMS];
02439       */
02440       int32 *dimsize = NULL;
02441       size_t var_type_size;
02442       int a;
02443         
02444       /* Add a variable to the end of the group's var list. */
02445       if ((retval = nc4_var_list_add(&grp->var, &var)))
02446         return retval;
02447       
02448       var->varid = grp->nvars++;
02449       var->created = NC_TRUE;
02450       var->written_to = NC_TRUE;
02451             
02452       /* Open this dataset in HDF4 file. */
02453       if ((var->sdsid = SDselect(h5->sdid, v)) == FAIL)
02454         return NC_EVARMETA;
02455 
02456       /* Get shape, name, type, and attribute info about this dataset. */
02457       if (!(var->name = malloc(NC_MAX_HDF4_NAME + 1)))
02458         return NC_ENOMEM;
02459       
02460       /* Invoke SDgetInfo with null dimsize to get rank. */
02461       if (SDgetinfo(var->sdsid, var->name, &rank, NULL, &data_type, &num_atts))
02462         return NC_EVARMETA;
02463       
02464       if(!(dimsize = (int32*)malloc(sizeof(int32)*rank)))
02465         return NC_ENOMEM;
02466       
02467       if (SDgetinfo(var->sdsid, var->name, &rank, dimsize, &data_type, &num_atts)) {
02468         if(dimsize) free(dimsize);
02469         return NC_EVARMETA;
02470       }
02471       
02472       var->ndims = rank;
02473       var->hdf4_data_type = data_type;
02474 
02475       /* Fill special type_info struct for variable type information. */
02476       if (!(var->type_info = calloc(1, sizeof(NC_TYPE_INFO_T)))) {
02477         if(dimsize) free(dimsize);
02478         return NC_ENOMEM;
02479       }
02480       
02481       if ((retval = get_netcdf_type_from_hdf4(h5, data_type, &var->type_info->nc_typeid, var->type_info))) {
02482         if(dimsize) free(dimsize);
02483         return retval;
02484       }
02485 
02486       /* Indicate that the variable has a pointer to the type */
02487       var->type_info->rc++;
02488       
02489       if ((retval = nc4_get_typelen_mem(h5, var->type_info->nc_typeid, 0, &var_type_size))) {
02490         if(dimsize) free(dimsize);
02491         return retval;
02492       }
02493 
02494       var->type_info->size = var_type_size;
02495       LOG((3, "reading HDF4 dataset %s, rank %d netCDF type %d", var->name, 
02496            rank, var->type_info->nc_typeid));
02497 
02498       /* Get the fill value. */
02499       if (!(var->fill_value = malloc(var_type_size))) {
02500         if(dimsize) free(dimsize);
02501         return NC_ENOMEM;
02502       }
02503 
02504       if (SDgetfillvalue(var->sdsid, var->fill_value))
02505       {
02506          /* Whoops! No fill value! */
02507          free(var->fill_value);
02508          var->fill_value = NULL;
02509       }
02510 
02511       /* Allocate storage for dimension info in this variable. */
02512       if (var->ndims)
02513       {
02514         if (!(var->dim = malloc(sizeof(NC_DIM_INFO_T *) * var->ndims))) {
02515           if(dimsize) free(dimsize);
02516           return NC_ENOMEM;
02517         }
02518         
02519         if (!(var->dimids = malloc(sizeof(int) * var->ndims))) {
02520           if(dimsize) free(dimsize);
02521           return NC_ENOMEM;
02522         }
02523       }
02524       
02525 
02526       /* Find its dimensions. */
02527       for (d = 0; d < var->ndims; d++)
02528       {
02529          int32 dimid, dim_len, dim_data_type, dim_num_attrs;
02530          char dim_name[NC_MAX_NAME + 1];
02531          NC_DIM_INFO_T *dim;
02532 
02533          if ((dimid = SDgetdimid(var->sdsid, d)) == FAIL) {
02534            if(dimsize) free(dimsize);
02535            return NC_EDIMMETA;
02536          }
02537          if (SDdiminfo(dimid, dim_name, &dim_len, &dim_data_type, 
02538                        &dim_num_attrs))
02539            {
02540              if(dimsize) free(dimsize);
02541              return NC_EDIMMETA;
02542            }
02543 
02544          /* Do we already have this dimension? HDF4 explicitly uses
02545           * the name to tell. */
02546          for (dim = grp->dim; dim; dim = dim->l.next)
02547             if (!strcmp(dim->name, dim_name))
02548                break;
02549 
02550          /* If we didn't find this dimension, add one. */
02551          if (!dim)
02552          {
02553             LOG((4, "adding dimension %s for HDF4 dataset %s", 
02554                  dim_name, var->name));
02555             if ((retval = nc4_dim_list_add(&grp->dim, &dim)))
02556                return retval;
02557             grp->ndims++;
02558             dim->dimid = grp->nc4_info->next_dimid++;
02559             if (strlen(dim_name) > NC_MAX_HDF4_NAME)
02560                return NC_EMAXNAME;
02561             if (!(dim->name = strdup(dim_name)))
02562                return NC_ENOMEM;
02563             if (dim_len)
02564                dim->len = dim_len;
02565             else
02566                dim->len = *dimsize;
02567          }
02568 
02569          /* Tell the variable the id of this dimension. */
02570          var->dimids[d] = dim->dimid;
02571       }
02572 
02573       /* Read the atts. */
02574       for (a = 0; a < num_atts; a++)
02575       {
02576          int32 att_data_type, att_count;
02577          size_t att_type_size;
02578 
02579          /* Add to the end of the list of atts for this var. */
02580          if ((retval = nc4_att_list_add(&var->att, &att))) {
02581            if(dimsize) free(dimsize);
02582            return retval;
02583          }
02584          att->attnum = var->natts++;
02585          att->created = NC_TRUE;
02586 
02587          /* Learn about this attribute. */
02588          if (!(att->name = malloc(NC_MAX_HDF4_NAME * sizeof(char)))) {
02589            if(dimsize) free(dimsize);
02590            return NC_ENOMEM;
02591          }
02592          if (SDattrinfo(var->sdsid, a, att->name, &att_data_type, &att_count)) {
02593            if(dimsize) free(dimsize);
02594             return NC_EATTMETA;
02595          }
02596          if ((retval = get_netcdf_type_from_hdf4(h5, att_data_type, 
02597                                                  &att->nc_typeid, NULL))) {
02598            if(dimsize) free(dimsize);
02599            return retval;
02600          }
02601          
02602          att->len = att_count;
02603 
02604          /* Allocate memory to hold the data. */
02605          if ((retval = nc4_get_typelen_mem(h5, att->nc_typeid, 0, &att_type_size))) {
02606            if(dimsize) free(dimsize);
02607            return retval;
02608          }
02609          if (!(att->data = malloc(att_type_size * att->len))) {
02610            if(dimsize) free(dimsize);
02611            return NC_ENOMEM;
02612          }
02613 
02614          /* Read the data. */
02615          if (SDreadattr(var->sdsid, a, att->data)) {
02616            if(dimsize) free(dimsize);
02617            return NC_EHDFERR;
02618          }
02619       }
02620       if(dimsize) free(dimsize);
02621    } /* next var */
02622 
02623 #ifdef LOGGING
02624    /* This will print out the names, types, lens, etc of the vars and
02625       atts in the file, if the logging level is 2 or greater. */ 
02626    log_metadata_nc(h5->root_grp->nc4_info->controller);
02627 #endif
02628    return NC_NOERR;   
02629    return NC_ENOTBUILT;
02630 }
02631 #endif /* USE_HDF4 */
02632 
02633 int
02634 NC4_open(const char *path, int mode, int basepe, size_t *chunksizehintp, 
02635          int use_parallel, void *mpidata, NC_Dispatch *dispatch, NC *nc_file)
02636 {
02637    int hdf_file = 0;
02638    MPI_Comm comm = MPI_COMM_WORLD;
02639    MPI_Info info = MPI_INFO_NULL;
02640    int res;
02641 
02642    assert(nc_file && path);
02643 
02644    LOG((1, "%s: path %s mode %d comm %d info %d", 
02645         __func__, path, mode, comm, info));
02646 
02647 #ifdef USE_PARALLEL
02648    if (mpidata) 
02649    { 
02650       comm = ((NC_MPI_INFO *)mpidata)->comm;
02651       info = ((NC_MPI_INFO *)mpidata)->info; 
02652    }
02653 #endif /* USE_PARALLEL */
02654     
02655    /* If this is our first file, turn off HDF5 error messages. */
02656    if (virgin)
02657    {
02658       if (H5Eset_auto(NULL, NULL) < 0)
02659          LOG((0, "Couldn't turn off HDF5 error messages!"));
02660       LOG((1, "HDF5 error messages turned off!"));
02661       virgin = 0;
02662    }
02663 
02664    /* Check the mode for validity. First make sure only certain bits
02665     * are turned on. Also MPI I/O and MPI POSIX cannot both be
02666     * selected at once. */
02667    if (mode & ~(NC_WRITE | NC_SHARE | NC_MPIIO | NC_MPIPOSIX | 
02668                 NC_PNETCDF | NC_NOCLOBBER | NC_NETCDF4 | NC_CLASSIC_MODEL) ||
02669        (mode & NC_MPIIO && mode & NC_MPIPOSIX))
02670       return NC_EINVAL;
02671 
02672 #ifndef USE_PARALLEL_POSIX
02673 /* If the HDF5 library has been compiled without the MPI-POSIX VFD, alias
02674  *      the NC_MPIPOSIX flag to NC_MPIIO. -QAK
02675  */
02676    if(mode & NC_MPIPOSIX)
02677    {
02678       mode &= ~NC_MPIPOSIX;
02679       mode |= NC_MPIIO;
02680    }
02681 #endif /* USE_PARALLEL_POSIX */
02682 
02683 
02684    /* Depending on the type of file, open it. */
02685 
02686 #if 0 /*def USE_PNETCDF*/
02687    if(mode & NC_PNETCDF) {
02688         /* this is not really an hdf file */
02689       int pnetcdf_nvars, i;
02690       NC_HDF5_FILE_INFO_T* nc4_info;
02691 
02692       /* Create the fake nc4_info data */
02693       res = nc4_nc4f_list_add(nc_file, path, mode);
02694 
02695       nc4_info = NC4_DATA(nc_file);
02696       assert(nc4_info);
02697 
02698       res = ncmpi_open(comm, path, mode, info, &(nc_file->int_ncid));
02699       nc4_info->pnetcdf_file++;
02700 
02701       /* Default to independent access, like netCDF-4/HDF5 files. */
02702       if (!res)
02703          res = ncmpi_begin_indep_data(nc_file->int_ncid);
02704 
02705       /* I need to keep track of the ndims of each var to translate
02706        * start, count, and stride arrays to MPI_Offset type. */
02707       if (!res)
02708       {
02709          res = ncmpi_inq_nvars(nc_file->int_ncid, &pnetcdf_nvars);
02710          for (i = 0; i < pnetcdf_nvars; i++)
02711             res = ncmpi_inq_varndims(nc_file->int_ncid, i, 
02712                                      &(nc4_info->pnetcdf_ndims[i]));
02713 
02714       }
02715    } else
02716 #endif
02717    {
02718       /* Figure out if this is a hdf4 or hdf5 file. */
02719      if ((res = nc_check_for_hdf(path, use_parallel, comm, info, &hdf_file)))
02720          return res;
02721 
02722       if (hdf_file == NC_HDF5_FILE)
02723       {
02724          nc_file->int_ncid = nc_file->ext_ncid;
02725          res = nc4_open_file(path, mode, comm, info, nc_file);
02726       }
02727 #ifdef USE_HDF4   
02728       else if (hdf_file == NC_HDF4_FILE)
02729       {
02730          nc_file->int_ncid = nc_file->ext_ncid;
02731          res = nc4_open_hdf4_file(path, mode, nc_file);
02732       }
02733 #endif /* USE_HDF4 */
02734       else /* netcdf */
02735       {
02736          assert(0); /* should never happen */
02737       }
02738    }
02739 
02740    return res;
02741 }
02742 
02743 /* Unfortunately HDF only allows specification of fill value only when
02744    a dataset is created. Whereas in netcdf, you first create the
02745    variable and then (optionally) specify the fill value. To
02746    accomplish this in HDF5 I have to delete the dataset, and recreate
02747    it, with the fill value specified. */
02748 /* QAK: This looks completely unused in the code. (?) */
02749 int 
02750 NC4_set_fill(int ncid, int fillmode, int *old_modep)
02751 {
02752    NC *nc;
02753    NC_HDF5_FILE_INFO_T* nc4_info;
02754  
02755    LOG((2, "%s: ncid 0x%x fillmode %d", __func__, ncid, fillmode));
02756 
02757    if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
02758       return NC_EBADID;
02759    assert(nc4_info);
02760 
02761    /* Trying to set fill on a read-only file? You sicken me! */
02762    if (nc4_info->no_write)
02763       return NC_EPERM;
02764 
02765    /* Did you pass me some weird fillmode? */
02766    if (fillmode != NC_FILL && fillmode != NC_NOFILL)
02767       return NC_EINVAL;
02768 
02769    /* If the user wants to know, tell him what the old mode was. */
02770    if (old_modep)
02771       *old_modep = nc4_info->fill_mode;
02772 
02773    nc4_info->fill_mode = fillmode;      
02774 
02775 #if 0 /*def USE_PNETCDF*/
02776    /* Take care of files created/opened with parallel-netcdf library. */
02777    if (nc4_info->pnetcdf_file)
02778      return ncmpi_set_fill(nc->int_ncid, fillmode, old_modep);
02779 #endif /* USE_PNETCDF */
02780 
02781 
02782    return NC_NOERR;
02783 }
02784 
02785 /* Put the file back in redef mode. This is done automatically for
02786  * netcdf-4 files, if the user forgets. */
02787 int
02788 NC4_redef(int ncid)
02789 {
02790   //NC *nc;
02791    NC_HDF5_FILE_INFO_T* nc4_info;
02792 
02793    LOG((1, "%s: ncid 0x%x", __func__, ncid));
02794 
02795    /* Find this file's metadata. */
02796    if (!(nc4_find_nc_file(ncid,&nc4_info)))
02797       return NC_EBADID;
02798    assert(nc4_info);
02799 
02800 #if 0 /*def USE_PNETCDF*/
02801    /* Take care of files created/opened with parallel-netcdf library. */
02802    if (nc4_info->pnetcdf_file)
02803       return ncmpi_redef(nc->int_ncid);
02804 #endif /* USE_PNETCDF */
02805 
02806    /* If we're already in define mode, return an error. */
02807    if (nc4_info->flags & NC_INDEF)
02808       return NC_EINDEFINE;
02809 
02810    /* If the file is read-only, return an error. */
02811    if (nc4_info->no_write)
02812       return NC_EPERM;
02813 
02814    /* Set define mode. */
02815    nc4_info->flags |= NC_INDEF;
02816 
02817    /* For nc_abort, we need to remember if we're in define mode as a
02818       redef. */
02819    nc4_info->redef = NC_TRUE;
02820 
02821    return NC_NOERR;
02822 }
02823 
02824 /* For netcdf-4 files, this just calls nc_enddef, ignoring the extra
02825  * parameters. */
02826 int
02827 NC4__enddef(int ncid, size_t h_minfree, size_t v_align,
02828             size_t v_minfree, size_t r_align)
02829 {
02830    if (nc4_find_nc_file(ncid,NULL) == NULL)
02831       return NC_EBADID;
02832 
02833    return NC4_enddef(ncid);
02834 }
02835 
02836 /* Take the file out of define mode. This is called automatically for
02837  * netcdf-4 files, if the user forgets. */
02838 static int NC4_enddef(int ncid)
02839 {
02840   NC *nc;
02841    NC_HDF5_FILE_INFO_T* nc4_info;
02842 
02843    LOG((1, "%s: ncid 0x%x", __func__, ncid));
02844    
02845    if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
02846       return NC_EBADID;
02847    assert(nc4_info);
02848 
02849 #if 0 /*def USE_PNETCDF*/
02850    if (nc4_info->pnetcdf_file)
02851    {
02852       int res;
02853       res = ncmpi_enddef(nc->int_ncid);
02854       if (!res)
02855       {
02856          if (nc4_info->pnetcdf_access_mode == NC_INDEPENDENT)
02857             res = ncmpi_begin_indep_data(nc->int_ncid);
02858       }
02859       return res;
02860    }
02861 #endif /* USE_PNETCDF */
02862 
02863    return nc4_enddef_netcdf4_file(nc4_info);
02864 }
02865 
02866 /* This function will write all changed metadata, and (someday) reread
02867  * all metadata from the file. */
02868 static int
02869 sync_netcdf4_file(NC_HDF5_FILE_INFO_T *h5)
02870 {
02871    int retval;
02872 
02873    assert(h5);
02874    LOG((3, "%s", __func__));
02875 
02876    /* If we're in define mode, that's an error, for strict nc3 rules,
02877     * otherwise, end define mode. */
02878    if (h5->flags & NC_INDEF)
02879    {
02880       if (h5->cmode & NC_CLASSIC_MODEL)
02881          return NC_EINDEFINE;
02882 
02883       /* Turn define mode off. */
02884       h5->flags ^= NC_INDEF;
02885       
02886       /* Redef mode needs to be tracked seperately for nc_abort. */
02887       h5->redef = NC_FALSE;
02888    }
02889 
02890 #ifdef LOGGING
02891    /* This will print out the names, types, lens, etc of the vars and
02892       atts in the file, if the logging level is 2 or greater. */ 
02893    log_metadata_nc(h5->root_grp->nc4_info->controller);
02894 #endif
02895 
02896    /* Write any metadata that has changed. */
02897    if (!(h5->cmode & NC_NOWRITE))
02898    {
02899       int bad_coord_order = 0;  /* if detected, propagate to all groups to consistently store dimids */
02900 
02901       if ((retval = nc4_rec_write_groups_types(h5->root_grp)))
02902          return retval;
02903       if ((retval = nc4_rec_detect_need_to_preserve_dimids(h5->root_grp, &bad_coord_order)))
02904          return retval;
02905       if ((retval = nc4_rec_write_metadata(h5->root_grp, bad_coord_order)))
02906          return retval;
02907    }
02908 
02909    if (H5Fflush(h5->hdfid, H5F_SCOPE_GLOBAL) < 0)
02910       return NC_EHDFERR;
02911 
02912    return retval;
02913 }
02914 
02915 /* Flushes all buffers associated with the file, after writing all
02916    changed metadata. This may only be called in data mode. */
02917 int
02918 NC4_sync(int ncid)
02919 {
02920    NC *nc;
02921    int retval;
02922    NC_HDF5_FILE_INFO_T* nc4_info;
02923 
02924    LOG((2, "%s: ncid 0x%x", __func__, ncid));
02925 
02926    if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
02927       return NC_EBADID;
02928    assert(nc4_info);
02929 
02930 #if 0 /*def USE_PNETCDF*/
02931    /* Take care of files created/opened with parallel-netcdf library. */
02932    if (nc4_info->pnetcdf_file)
02933       return ncmpi_sync(nc->int_ncid);
02934 #endif /* USE_PNETCDF */
02935 
02936    /* If we're in define mode, we can't sync. */
02937    if (nc4_info && nc4_info->flags & NC_INDEF)
02938    {
02939       if (nc4_info->cmode & NC_CLASSIC_MODEL)
02940          return NC_EINDEFINE;
02941       if ((retval = NC4_enddef(ncid)))
02942          return retval;
02943    }
02944 
02945    return sync_netcdf4_file(nc4_info);
02946 }
02947 
02948 /* This function will free all allocated metadata memory, and close
02949    the HDF5 file. The group that is passed in must be the root group
02950    of the file. */
02951 static int
02952 close_netcdf4_file(NC_HDF5_FILE_INFO_T *h5, int abort)
02953 {
02954    int retval = NC_NOERR;
02955 
02956    assert(h5 && h5->root_grp);
02957    LOG((3, "%s: h5->path %s abort %d", __func__, h5->controller->path, abort));
02958 
02959    /* According to the docs, always end define mode on close. */
02960    if (h5->flags & NC_INDEF)
02961       h5->flags ^= NC_INDEF;
02962 
02963    /* Sync the file, unless we're aborting, or this is a read-only
02964     * file. */
02965    if (!h5->no_write && !abort)
02966       if ((retval = sync_netcdf4_file(h5)))
02967         goto exit;
02968 
02969    /* Delete all the list contents for vars, dims, and atts, in each
02970     * group. */
02971    if ((retval = nc4_rec_grp_del(&h5->root_grp, h5->root_grp)))
02972         goto exit;
02973 
02974    /* Close hdf file. */
02975 #ifdef USE_HDF4
02976    if (h5->hdf4)
02977    {
02978       if (SDend(h5->sdid))
02979          BAIL_QUIET(NC_EHDFERR);
02980    } 
02981    else
02982 #endif /* USE_HDF4 */
02983    {
02984 #ifdef USE_PARALLEL
02985       /* Free the MPI Comm & Info objects, if we opened the file in parallel */
02986       if(h5->parallel)
02987       {
02988           MPI_Comm_free(&h5->comm);
02989           if(MPI_INFO_NULL != h5->info)
02990               MPI_Info_free(&h5->info);
02991       }
02992 #endif
02993       if (H5Fclose(h5->hdfid) < 0) 
02994       {
02995         int nobjs;
02996 
02997         nobjs = H5Fget_obj_count(h5->hdfid, H5F_OBJ_ALL);
02998         /* Apparently we can get an error even when nobjs == 0 */
02999         if(nobjs < 0) {
03000           BAIL_QUIET(NC_EHDFERR);
03001         } else if(nobjs > 0) {
03002 #ifdef LOGGING
03003          /* If the close doesn't work, probably there are still some HDF5
03004           * objects open, which means there's a bug in the library. So
03005           * print out some info on to help the poor programmer figure it
03006           * out. */
03007          LOG((0, "There are %d HDF5 objects open!", nobjs));
03008 #endif      
03009          BAIL_QUIET(NC_EHDFERR);
03010         }
03011       }
03012    }
03013 
03014 exit:
03015    /* Free the nc4_info struct; above code should have reclaimed 
03016       everything else */
03017    if(h5 != NULL)
03018        free(h5);
03019 
03020    return retval;
03021 }
03022 
03023 /* From the netcdf-3 docs: The function nc_abort just closes the
03024    netCDF dataset, if not in define mode. If the dataset is being
03025    created and is still in define mode, the dataset is deleted. If
03026    define mode was entered by a call to nc_redef, the netCDF dataset
03027    is restored to its state before definition mode was entered and the
03028    dataset is closed. */
03029 int
03030 NC4_abort(int ncid)
03031 {
03032    NC *nc;
03033    int delete_file = 0;
03034    char path[NC_MAX_NAME + 1];
03035    int retval = NC_NOERR;
03036    NC_HDF5_FILE_INFO_T* nc4_info;
03037 
03038    LOG((2, "%s: ncid 0x%x", __func__, ncid));
03039 
03040    /* Find metadata for this file. */
03041    if (!(nc = nc4_find_nc_file(ncid,&nc4_info)))
03042       return NC_EBADID;
03043 
03044    assert(nc4_info);
03045 
03046 #if 0 /*def USE_PNETCDF*/
03047    /* Take care of files created/opened with parallel-netcdf library. */
03048    if (nc4_info->pnetcdf_file)
03049       return ncmpi_abort(nc->int_ncid);
03050 #endif /* USE_PNETCDF */
03051 
03052    /* If we're in define mode, but not redefing the file, delete it. */
03053    if (nc4_info->flags & NC_INDEF && !nc4_info->redef)
03054    {
03055       delete_file++;
03056       strncpy(path, nc->path,NC_MAX_NAME);
03057    }
03058 
03059    /* Free any resources the netcdf-4 library has for this file's
03060     * metadata. */
03061    if ((retval = close_netcdf4_file(nc4_info, 1)))
03062       return retval;
03063    
03064    /* Delete the file, if we should. */
03065    if (delete_file)
03066       if (remove(path) < 0)
03067           return NC_ECANTREMOVE;
03068 
03069    return retval;
03070 }
03071 
03072 /* Close the netcdf file, writing any changes first. */
03073 int
03074 NC4_close(int ncid)
03075 {
03076    NC_GRP_INFO_T *grp;
03077    NC *nc;
03078    NC_HDF5_FILE_INFO_T *h5;
03079    int retval;
03080 
03081    LOG((1, "%s: ncid 0x%x", __func__, ncid));
03082 
03083    /* Find our metadata for this file. */
03084    if ((retval = nc4_find_nc_grp_h5(ncid, &nc, &grp, &h5)))
03085       return retval;
03086 
03087    assert(nc && h5 && grp);
03088 
03089    /* This must be the root group. */
03090    if (grp->parent)
03091       return NC_EBADGRPID;
03092 
03093 #if 0 /*def USE_PNETCDF*/
03094    /* Take care of files created/opened with parallel-netcdf library. */
03095    if (h5->pnetcdf_file)
03096       return ncmpi_close(nc->int_ncid);
03097 #endif /* USE_PNETCDF */
03098 
03099    /* Call the nc4 close. */
03100    if ((retval = close_netcdf4_file(grp->nc4_info, 0)))
03101       return retval;
03102 
03103    return NC_NOERR;
03104 }
03105 
03106 /* It's possible for any of these pointers to be NULL, in which case
03107    don't try to figure out that value. */
03108 int
03109 NC4_inq(int ncid, int *ndimsp, int *nvarsp, int *nattsp, int *unlimdimidp)
03110 {
03111    NC *nc;
03112    NC_HDF5_FILE_INFO_T *h5;
03113    NC_GRP_INFO_T *grp;
03114    NC_DIM_INFO_T *dim;
03115    NC_ATT_INFO_T *att;
03116    NC_VAR_INFO_T *var;
03117    int retval;
03118 
03119    LOG((2, "%s: ncid 0x%x", __func__, ncid)); 
03120 
03121    /* Find file metadata. */
03122    if ((retval = nc4_find_nc_grp_h5(ncid, &nc, &grp, &h5)))
03123       return retval;
03124 
03125    assert(h5 && grp && nc);
03126 
03127 #if 0 /*def USE_PNETCDF*/
03128    /* Take care of files created/opened with parallel-netcdf library. */
03129    if (h5->pnetcdf_file)
03130       return ncmpi_inq(nc->int_ncid, ndimsp, nvarsp, nattsp, unlimdimidp);
03131 #endif /* USE_PNETCDF */
03132 
03133    /* Count the number of dims, vars, and global atts. */
03134    if (ndimsp)
03135    {
03136       *ndimsp = 0;
03137       for (dim = grp->dim; dim; dim = dim->l.next)
03138          (*ndimsp)++;
03139    }
03140    if (nvarsp)
03141    {
03142       *nvarsp = 0;
03143       for (var = grp->var; var; var= var->l.next)
03144         (*nvarsp)++;
03145    }
03146    if (nattsp)
03147      {
03148       *nattsp = 0;
03149       for (att = grp->att; att; att = att->l.next)
03150          (*nattsp)++;
03151    }
03152 
03153    if (unlimdimidp)
03154    {
03155       /* Default, no unlimited dimension */
03156       *unlimdimidp = -1;
03157 
03158       /* If there's more than one unlimited dim, which was not possible
03159          with netcdf-3, then only the last unlimited one will be reported
03160          back in xtendimp. */
03161       /* Note that this code is inconsistent with nc_inq_unlimid() */
03162       for (dim = grp->dim; dim; dim = dim->l.next)
03163          if (dim->unlimited)
03164          {
03165             *unlimdimidp = dim->dimid;
03166             break;
03167          }
03168    }
03169 
03170    return NC_NOERR;   
03171 }
03172 
03173 
03174 /* This function will do the enddef stuff for a netcdf-4 file. */
03175 int
03176 nc4_enddef_netcdf4_file(NC_HDF5_FILE_INFO_T *h5)
03177 {
03178    assert(h5);
03179    LOG((3, "%s", __func__));
03180 
03181    /* If we're not in define mode, return an error. */
03182    if (!(h5->flags & NC_INDEF))
03183       return NC_ENOTINDEFINE;
03184 
03185    /* Turn define mode off. */
03186    h5->flags ^= NC_INDEF;
03187 
03188    /* Redef mode needs to be tracked seperately for nc_abort. */
03189    h5->redef = NC_FALSE;
03190 
03191    return sync_netcdf4_file(h5);
03192 }
03193 
03194 #ifdef EXTRA_TESTS
03195 int
03196 nc_exit()
03197 {
03198    if (num_plists || num_spaces)
03199       return NC_EHDFERR;
03200       
03201    return NC_NOERR;
03202 }
03203 #endif /* EXTRA_TESTS */
03204 
03205 #ifdef USE_PARALLEL
03206 int
03207 nc_use_parallel_enabled()
03208 {
03209    return 0;
03210 }
03211 #endif /* USE_PARALLEL */
03212 
03213 
 All Data Structures Files Functions Variables Typedefs Defines