HDFEOS2.h

Go to the documentation of this file.

//  This file is part of the hdf4 data handler for the OPeNDAP data server.
//
// Author:   Kent Yang,Choonghwan Lee <myang6@hdfgroup.org>
// Copyright (c) 2010-2012 The HDF Group

//#include "InternalErr.h"
#ifdef USE_HDFEOS2_LIB

#ifndef _HDFEOS2_H
#define _HDFEOS2_H
#ifdef _WIN32
#ifdef _DEBUG
#define _CRTDBG_MAP_ALLOC
#include <stdlib.h>
#include <crtdbg.h>
#endif
#endif


#include <iostream>
#include <string>
#include <vector>
#include <map>
#include <set>
#include "mfhdf.h"
#include "hdf.h"
#include "HdfEosDef.h"


// Add MISR SOM projection header
#include "misrproj.h"

#include "HDFEOS2EnumType.h"

#ifdef _WIN32
#ifdef _MSC_VER
#pragma warning(disable:4290)
#endif
#endif

namespace HDFEOS2
{
    class Exception:public std::exception
    {
        public:
            Exception (const std::string & msg)
                : message (msg), isHDFEOS2 (true)
            {
            }

            virtual ~ Exception () throw ()
            {
            }

            virtual const char *what () const throw ()
            {
                return this->message.c_str ();
            }

            virtual bool getFileType ()
            {
                return this->isHDFEOS2;
            }

            virtual void setFileType (bool isHDFEOS2)
            {
                this->isHDFEOS2 = isHDFEOS2;
            }

            virtual void setException (std::string message)
            {
                this->message = message;
            }

        protected:
            std::string message;
            bool isHDFEOS2;
    };

    template < typename T > class LightVector {
        public:
            LightVector ()
                : data (0), length (0), capacity (0) {
            }

            LightVector (const LightVector < T > &that)
            {
                this->data = new T[that.length];
                for (unsigned int i = 0; i < that.length; ++i)
                    this->data[i] = that[i];
                this->length = that.length;
                this->capacity = that.length;
            }

            ~LightVector () {
                if (this->data)
                    delete[]data;
            }

            void push_back (const T & d)
            {
                this->reserve (this->length + 1);
                this->data[this->length] = d;
                ++this->length;
            }

            void reserve (unsigned int len)
            {
                if (this->capacity >= len)
                    return;

                this->capacity = len;
                T *old = this->data;

                this->data = new T[len];
                if (old) {
                    for (unsigned int i = 0; i < this->length; ++i)
                        this->data[i] = old[i];
                    delete[]old;
                }
            }

            void resize (unsigned int len)
            {
                if (this->length == len)
                    return;
                else if (this->length < len) {
                    if (this->capacity < len) {
                        this->capacity = len;
                        T *old = this->data;

                        this->data = new T[len];
                        if (old) {
                            for (unsigned int i = 0; i < this->length; ++i)
                                this->data[i] = old[i];
                            delete[]old;
                        }
                    }
                }
                else {
                    this->capacity = len;
                    T *old = this->data;

                    this->data = new T[len];
                    for (unsigned int i = 0; i < len; ++i)
                        this->data[i] = old[i];
                    if (old)
                        delete[]old;
                }
                this->length = len;
            }

            unsigned int size () const
            {
                return this->length;
            }

            T & operator[] (unsigned int i)
            {
                return this->data[i];
            }
            const T & operator[] (unsigned int i) const
            {
                return this->data[i];
            }

            LightVector < T > &operator= (const LightVector < T > &that)
            {
                if (this != &that) {
                    this->data = new T[that.length];
                    for (unsigned int i = 0; i < that.length; ++i)
                        this->data[i] = that[i];
                    this->length = that.length;
                    this->capacity = that.length;
                }
                return *this;
            }

        private:
            T * data;
            unsigned int length;
            unsigned int capacity;
    };

    class SwathDimensionAdjustment;

    class Dimension
    {
        public:
            const std::string & getName () const
            {
                return this->name;
            }
            int32 getSize () const
            {
                return this->dimsize;
            }

        protected:
            Dimension (const std::string & name, int32 dimsize)
                : name (name), dimsize (dimsize)
            {
            }

        protected:
            std::string name;
            int32 dimsize;

        friend class File;
        friend class Dataset;
        friend class SwathDimensionAdjustment;
    };

    class Field
    {
        public:
            Field ()
                :fieldtype (0), condenseddim (false), iscoard (false), ydimmajor (true), speciallon (false), specialformat (0), haveaddedfv (false), addedfv (-9999.0), dmap (false)/*, field_cache(0)*/

            {
            }
            virtual ~ Field ();

        public:

            const std::string & getName () const
            {
                return this->name;
            }

            const std::string & getNewName () const
            {
                return this->newname;
            }

            int32 getRank () const 
            {   
                return this->rank;
            }

            int32 getType () const
            {
                return this->type;
            }

            const std::vector < Dimension * >&getCorrectedDimensions () const
            {
                return this->correcteddims;
            }

            std::vector < Dimension * >*getCorrectedDimensionsPtr ()
            {
                return &(this->correcteddims);
            }

            void setCorrectedDimensions (std::vector < Dimension * >dims)
            {
                correcteddims = dims;
            }

            const std::string getCoordinate () const
            {
                return this->coordinates;
            }

            void setCoordinates (std::string coor)
            {
                coordinates = coor;
            }

            const std::string getUnits () const
            {
                return this->units;
            }

            void setUnits (std::string uni)
            {
                units = uni;
            }

            const float getAddedFillValue () const
            {
                return this->addedfv;
            }

            // Add the "_FillValue" attribute
            // This is for supporting the old versions of HDF-EOS2 data(AIRS) because
            // some data products have fillvalue(-9999.0) but don't specify the fillvalue.
            // We add the fillvalue to ensure the netCDF client can successfully display the data.
            // KY 2013-06-30
            void addFillValue (float fv)
            {
                addedfv = fv;
            }

            const bool haveAddedFillValue () const
            {
                return this->haveaddedfv;
            }

            void setAddedFillValue (bool havefv)
            {
                haveaddedfv = havefv;
            }

            
            // Obtain fieldtype values
            // For fieldtype values:
            // 0 is general fields
            // 1 is latitude.
            // 2 is longtitude.    
            // 3 is defined level.
            // 4 is an inserted natural number.
            // 5 is time.

            const int getFieldType () const
            {
                return this->fieldtype;
            }

            const std::vector < Dimension * >&getDimensions () const
            {
                return this->dims;
            }

            const std::vector < char >&getFillValue () const
            {
                return this->filler;
            }

            const bool getYDimMajor () const
            {
                return this->ydimmajor;
            }

            const bool getSpecialLon () const
            {
                return this->speciallon;
            }

            const int getSpecialLLFormat () const
            {
                return this->specialformat;
            }

            const bool getCondensedDim () const
            {
                return this->condenseddim;
            }

            const bool UseDimMap () const
            {
                return this->dmap;
            }

//No need to check field_cache in the DDS level. 
//May remove the debugging info. totally in the next release.
//KY 2014-10-23
#if 0
            // 2 exactly cached 
            // 1 maybe cached, need to check the 
            // file size when accessing the data
            const short UseFieldCache () const
            {
                return this->field_cache;
            }
#endif
        protected:

            // field name
            std::string name;

            // field dimension rank
            int32 rank;

            // field  datatype
            int32 type;

            // field dimensions with original dimension names
            std::vector < Dimension * >dims;

            // field dimensions with the corrected(CF) dimension names
            std::vector < Dimension * >correcteddims;
            
            // This is for reading the fillvalue. 
            // HDF-EOS2 provides a special routine to read fillvalue.
            // Up to HDF-EOS2 version 2.18, this is the only field attribute
            // that HDF-EOS2 APIs provide. KY 2013-07-01
            
            std::vector < char >filler;

            // Coordinate attributes that includes coordinate variable list.
            std::string coordinates;

            // newname is to record CF Grid/Swath name + "_"+ CF field name(special characters replaced by underscores).
            std::string newname;


            // This flag will specify the fieldtype.
            // 0 means this field is general field.
            // 1 means this field is lat.
            // 2 means this field is lon.
            // 3 means this field is other dimension variable.
            // 4 means this field is added other dimension variable with nature number.
            // 5 means time, but currently the units is not correct.
            int fieldtype;

            //  Latitude and longitude retrieved by HDF-EOS2 are always
            //  2-D arrays(XDim * YDim). However, for some projections
            // (geographic etc.), latiude and longitude can be condensed to
            //  1-D arrays. The handler will track such projections and condense
            //  the latitude and longitude to 1-D arrays. This can reduce       
            //  the disk storage and can greatly improve the performance for 
            //  the visualization tool to access the latitude and longitde.
            //  condenseddim is the flag internally used by the handler to track this.
            bool condenseddim;

            //   This flag is to mark if the data should follow COARDS.
            bool iscoard;

            //   This flag is to check if the field is YDim major(temp(YDim,XDim). This
            //   flag is necessary when calling GDij2ll to retrieve latitude and longitude.
            bool ydimmajor;

            // SOme special longitude is from 0 to 360.We need to check this case with this flag.
            bool speciallon;

            // This flag specifies the special latitude/longitude coordinate format
            // The latiude and longitude should represent as DDDMMMSSS format
            // However, we found some files simply represent lat/lon as -180.0000000 or -90.000000. 
            // Some other files use default. So we this flag to record this and correctly retrieve
            // the latitude and longitude values in the DAP output.
            // 0 means normal 
            // 1 means the coordinate is -180 to 180
            // 2 means the coordinate is default(0)
            int specialformat;

            // CF units attribute(mostly to add latitude and longitude CF units).
            std::string units;

            // Some data products have fillvalue(-9999.0) but don't specify the fillvalue.
            // We add the fillvalue to ensure the netCDF client can successfully display the data.
            // haveaddedfv and addedfv are to check if having added fillvalues.
            bool haveaddedfv;
            float addedfv;

            // Check if this swath uses the dimension map. 
            bool dmap;

            // Flag to indicate the field cache. Now this only applies to some HDF-EOS2 grid lat/lon fields.
            //short field_cache;

            // This may not necessary. Check in the next release. KY 2012-09-20
            // bool need_adjustment;

        friend class Dataset;
        friend class SwathDimensionAdjustment;
        friend class GridDataset;
        friend class SwathDataset;
        friend class File;
    };

#if 0
    // For future improvement of the modulization
    class GeoField:public Field
    {

        protected:
            bool condenseddim;
            bool ydimmajor;
            bool speciallon;

    };
#endif

    class Attribute
    {
        public:

            const std::string & getName () const
            {
                return this->name;
            }

            const std::string & getNewName () const
            {
                return this->newname;
            }

            int32 getType () const
            {
                return this->type;
            }

            int32 getCount () const
            {
                return this->count;
            }

            const std::vector < char >&getValue () const
            {
                return this->value;
            }

        protected:

            std::string name;

            std::string newname;

            int32 type;

            int32 count;

            std::vector < char >value;

        friend class Dataset;
    };

    class Dataset
    {
        public:
            const std::string & getName () const
            {
                return this->name;
            }
            const std::vector < Dimension * >&getDimensions () const
            {
                return this->dims;
            }
            const std::vector < Field * >&getDataFields () const
            {
                return this->datafields;
            }

            const std::vector < Attribute * >&getAttributes () const
            {
                return this->attrs;
            }

            const SOType getScaleType () const
            {
                return this->scaletype;
            }


        protected:
            Dataset (const std::string & n)
                : datasetid (-1), addfvalueattr(false),name (n),scaletype(DEFAULT_CF_EQU)
            {
            }

            virtual ~ Dataset ();

            void ReadDimensions (int32 (*entries) (int32, int32, int32 *),
                int32 (*inq) (int32, char *, int32 *),
                std::vector < Dimension * >&dims) throw (Exception);

            void ReadFields (int32 (*entries) (int32, int32, int32 *),
                int32 (*inq) (int32, char *, int32 *, int32 *),
                intn (*fldinfo) (int32, char *, int32 *, int32 *,
                int32 *, char *),
                intn (*readfld) (int32, char *, int32 *, int32 *,
                int32 *, VOIDP),
                intn (*getfill) (int32, char *, VOIDP),
                bool geofield, std::vector < Field * >&fields) 
                throw (Exception);

            void ReadAttributes (int32 (*inq) (int32, char *, int32 *),
                intn (*attrinfo) (int32, char *, int32 *, int32 *),
                intn (*readattr) (int32, char *, VOIDP),
                std::vector < Attribute * >&attrs)
                throw (Exception);
 
            void SetScaleType(const std::string EOS2ObjName) throw(Exception);

        protected:
            int32 datasetid;

            bool addfvalueattr;

            std::string name;

            std::vector < Dimension * >dims;

            std::vector < Field * >datafields;

            std::vector < Attribute * >attrs;

            std::map < std::string, std::string > dimcvarlist;

            std::map < std::string, std::string > ncvarnamelist;

            std::map < std::string, std::string > ndimnamelist;

            // Some MODIS files don't use the CF linear equation y = scale * x + offset,
            // The scaletype distinguishs products following different scale and offset rules.
            // Note the assumption here: we assume that all fields will
            // use one scale and offset function in a file. If
            // multiple scale and offset equations are used in one file, our
            // function will fail. So far only one scale and offset equation is
            //  applied for NASA HDF-EOS2 files we observed. KY 2012-6-13
            // Since I found one MODIS product(MOD09GA) uses different scale offset
            // equations for different grids. I had to move the scaletype to the
            // group level. Hopefully this is the final fix. Truly hope that 
            // this will not happen at the field level since it will be too messy to 
            // check.  KY 2012-11-21
            SOType scaletype;

        friend class File;
    };

    class GridDataset:public Dataset
    {
        public:
            class Info
            {
                public:

                    int32 getX ()const
                    {
                        return this->xdim;
                    }

                    int32 getY () const
                    {
                        return this->ydim;
                    }

                    const float64 *getUpLeft () const
                    {
                        return this->upleft;
                    }

                    const float64 *getLowRight () const
                    {
                        return this->lowright;
                    }

                protected:
                    int32 xdim;
                    int32 ydim;
                    float64 upleft[2];
                    float64 lowright[2];

                friend class GridDataset;
            };

            class Projection
            {
                public:

                    int32 getCode () const
                    {
                        return this->code;
                    }

                    int32 getZone () const
                    {
                        return this->zone;
                    }

                    int32 getSphere () const
                    {
                        return this->sphere;
                    }

                    const float64 *getParam () const
                    {
                        return this->param;
                    }

                    int32 getPix () const
                    {
                        return this->pix;
                    }

                    int32 getOrigin () const
                    {
                        return this->origin;
                    }

                protected:
                    int32 code;
                    int32 zone;
                    int32 sphere;
                    float64 param[16];
                    int32 pix;
                    int32 origin;

                friend class GridDataset;
            };

            class Calculated
            {
                public:

                    bool isYDimMajor () throw (Exception);

                protected:

                    Calculated (const GridDataset * grid)
                        : grid (grid),  ydimmajor (false)
                    {
                    }

                    Calculated & operator= (const Calculated & victim)
                    {
                        if (this != &victim) {
                            this->grid = victim.grid;
                            this->ydimmajor = victim.ydimmajor;
                        }
                        return *this;
                    }

                    void DetectMajorDimension () throw (Exception);

                    int DetectFieldMajorDimension () throw (Exception);


                protected:
                    const GridDataset *grid;
                    bool ydimmajor;

                friend class GridDataset;
                friend class File;
            };

            public:
                static GridDataset *Read (int32 fd, const std::string & gridname) throw (Exception);

                virtual ~ GridDataset ();

                const Info & getInfo () const
                {
                    return this->info;
                }

                const Projection & getProjection () const
                {
                    return this->proj;
                }

                Calculated & getCalculated () const;

                void setDimxName (std::string dxname)
                {
                    dimxname = dxname;
                }
                void setDimyName (std::string dyname)
                {
                    dimyname = dyname;
                }

                const bool getLatLonFlag () const
                {
                    return this->ownllflag;
                }

            private:
                GridDataset (const std::string & name)
                    : Dataset (name), calculated (0), ownllflag (false), iscoard (false)
                {
                }

            protected:

                Info info;

                Projection proj;

                mutable Calculated calculated;

                bool ownllflag;

                bool iscoard;

                Field *latfield;
                Field *lonfield;

                std::string dimxname;
                std::string dimyname;

            friend class File;

        };

        class File;


        class SwathDataset:public Dataset
        {

            public:

                class DimensionMap
                {
                    public:
                        const std::string & getGeoDimension () const
                        {
                            return this->geodim;
                        }
                        const std::string & getDataDimension () const
                        {
                            return this->datadim;
                        }
                        int32 getOffset () const
                        {
                            return this->offset;
                        }
                        int32 getIncrement () const
                        {
                            return this->increment;
                        }

                    protected:
                        DimensionMap (const std::string & geodim, const std::string & datadim, int32 offset, int32 increment)
                            : geodim (geodim), datadim (datadim), offset (offset), increment (increment)
                        {
                        }

                    protected:

                        std::string geodim;
                        std::string datadim;
                        int32 offset;
                        int32 increment;

                    friend class SwathDataset;
                    friend class SwathDimensionAdjustment;
                    friend class File;
                };

                class IndexMap
                {
                    public:
                        const std::string & getGeoDimension () const
                        {
                            return this->geo;
                        }
                        const std::string & getDataDimension () const
                        {
                            return this->data;
                        }
                        const LightVector < int32 > &getIndices () const
                        {
                            return this->indices;
                        }

                    protected:
                        std::string geo;
                        std::string data;
                        LightVector < int32 > indices;

                    friend class SwathDataset;
                };

            public:
                static SwathDataset *Read (int32 fd, const std::string & swathname) throw (Exception);

                virtual ~ SwathDataset ();

                const std::vector < DimensionMap * >&getDimensionMaps () const
                {
                    return this->dimmaps;
                }
                const std::vector < IndexMap * >&getIndexMaps () const
                {
                    return this->indexmaps;
                }

                const std::vector < Field * >&getGeoFields () const
                {
                    return this->geofields;
                }


                void set_num_map (int this_num_map)
                {
                    num_map = this_num_map;
                }
                int get_num_map () const
                {
                    return num_map;
                };

            private:
                 SwathDataset (const std::string & name)
                    : Dataset (name) {
                 }


                int ReadDimensionMaps (std::vector < DimensionMap * >&dimmaps) throw (Exception);

                void ReadIndexMaps (std::vector < IndexMap * >&indexmaps) throw (Exception);

            protected:

                std::vector < DimensionMap * >dimmaps;

                std::vector < IndexMap * >indexmaps;

                std::set < std::string > nonmisscvdimlist;

                std::vector < Field * >geofields;

                int num_map;

            friend class File;
        };

        class PointDataset:public Dataset
        {
            public:
                static PointDataset *Read (int32 fd, const std::string & pointname) throw (Exception);
                virtual ~ PointDataset ();

            private:
                PointDataset (const std::string & name)
                    : Dataset (name)
                {
                }
        };


        class File
        {
            public:

                static File *Read (const char *path,int32 gridfd,int32 swathfd) throw (Exception);


                void Prepare(const char *path) throw(Exception);

                bool check_special_1d_grid() throw(Exception);


                bool getOneLatLon ()
                {
                    return this->onelatlon;
                }

                ~File ();

                const std::string & getPath () const
                {
                    return this->path;
                }

                const std::vector < GridDataset * >&getGrids () const
                {
                    return this->grids;
                }

                const std::vector < SwathDataset * >&getSwaths () const
                {
                    return this->swaths;
                }

                const std::vector < PointDataset * >&getPoints () const
                {
                    return this->points;
                }

                const std::string get_first_grid_name() const 
                { return this->grids[0]->getName();}


            protected:
                File (const char *path)
                    : path (path), onelatlon (false), iscoard (false), gridfd (-1), swathfd (-1)
                {
                }

            protected:

                std::string path;

                std::vector < GridDataset * >grids;

                std::vector < SwathDataset * >swaths;

                std::vector < PointDataset * >points;

                // This is for the case that only one lat/lon
                // set is provided for multiple grid cases.
                //  The current case is that the user provides
                // the latitude and longitude under a special grid.
                // By default, the grid name is "location". This will
                // cover the AIRS grid case.
                bool onelatlon;

                bool iscoard;

                /* 
                  A grid's X-dimension can have different names: XDim, LatDim, etc.
                 * Y-dimension also has YDim, LonDim, etc.
                 * This function returns the name of X-dimension which is used in
                 * the given file.
                 * For better performance, we check the first grid or swath only.
                 */
                std::string get_geodim_x_name ();
                std::string get_geodim_y_name ();

                // Internal function used by  
                // get_geodim_x_name and get_geodim_y_name functions.
                // This function is not intended to be used outside the 
                // get_geodim_x_name and get_geodim_y_name functions.
                void _find_geodim_names ();

                std::string _geodim_x_name;
                std::string _geodim_y_name;
                static const char *_geodim_x_names[];
                static const char *_geodim_y_names[];

                std::string get_latfield_name ();
                std::string get_lonfield_name ();

                // Internal function used by  
                // get_latfield_name and get_lonfield_name functions.
                // This function is not intended to be used outside 
                // the get_latfield_name and get_lonfield_name functions.
                void _find_latlonfield_names ();

                std::string _latfield_name;
                std::string _lonfield_name;
                static const char *_latfield_names[];
                static const char *_lonfield_names[];

                std::string get_geogrid_name ();

                // Internal function used by
                // the get_geogrid_name function.
                // This function is not intended to be used outside the get_geogrid_name function.
                void _find_geogrid_name ();

                std::string _geogrid_name;
                static const char *_geogrid_names[];


                // All the following functions are called by the Prepare() function.

                // Check if we have the dedicated lat/lon grid.
                void check_onelatlon_grids();                

                // For one grid, need to handle the third-dimension(both existing and missing) coordinate variables
                void handle_one_grid_zdim(GridDataset*);

                // For one grid, need to handle lat/lon(both existing lat/lon and calculated lat/lon from EOS2 APIs)
                void handle_one_grid_latlon(GridDataset*) throw(Exception);

                // For the case of which all grids have one dedicated lat/lon grid,
                // this function shows how to handle lat/lon fields.
                void handle_onelatlon_grids() throw (Exception);

                // Handle the dimension name to coordinate variable map for grid. 
                void handle_grid_dim_cvar_maps() throw(Exception);

                // Follow COARDS for grids.
                void handle_grid_coards() throw(Exception);

                // Create the corrected dimension vector for each field when COARDS is not followed.
                void update_grid_field_corrected_dims() throw(Exception);

                // Handle CF attributes for grids. 
                // The CF attributes include "coordinates", "units" for coordinate variables and "_FillValue". 
                void handle_grid_cf_attrs() throw(Exception);

                // Special handling SOM(Space Oblique Mercator) projection files
                void handle_grid_SOM_projection() throw(Exception);

                // Obtain the number of dimension maps in this file. The input parameter is the number of swath.
                int  obtain_dimmap_num(int numswath) throw(Exception);

                // Create the dimension name to coordinate variable name map for lat/lon. 
                // The input parameter is the number of dimension maps in this file.
                void create_swath_latlon_dim_cvar_map(int numdm) throw(Exception);

                // Create the dimension name to coordinate variable name map for non lat/lon coordinate variables.
                void create_swath_nonll_dim_cvar_map() throw(Exception);

                // Handle swath dimension name to coordinate variable name maps. 
                // The input parameter is the number of dimension maps in this file.
                void handle_swath_dim_cvar_maps(int numdm) throw(Exception);

                // Handle CF attributes for swaths. 
                // The CF attributes include "coordinates", "units" for coordinate variables and "_FillValue". 
                void handle_swath_cf_attrs() throw(Exception);

             private:

                // HDF-EOS2 Grid File ID. Notice this ID is not an individual grid ID but the grid file ID returned by 
                // calling the HDF-EOS2 API GDopen. 
                int32 gridfd;

                // HDF-EOS2 Swath File ID. Notice this ID is not an individual swath ID but the swath file ID returned by
                 // calling the HDF-EOS2 API SWopen.
                int32 swathfd;


        };


        struct Utility
        {

            static bool ReadNamelist (const char *path,
                int32 (*inq) (char *, char *, int32 *),
                std::vector < std::string > &names);


        };

}
#endif

#endif