source: XMLIO_V2/dev/dev_rv/src/XMLIO/NetCDF4_data_output.hpp @ 127

#ifndef __NETCDF4_DATA_OUTPUT__
#define __NETCDF4_DATA_OUTPUT__

#include <netcdfcpp.h>

namespace XMLIOSERVER
{

   static const char* TimeName = "time";

   class NetCDF4DataOutput : public AbstractDataOutput
   {
      public :

         NetCDF4DataOutput(CFile* const _file)
            : AbstractDataOutput(_file), dataFile(NULL)
         { /* Ne rien faire de plus */ }

         const NcFile* getDataFile(void) const { return (dataFile); }

         void writeVarData(const string& id, const Array<float, 1>& vdata)
         {
            //if (!latVar->put(lats, NLAT))
            //   return NC_ERR;
         }

         virtual ~NetCDF4DataOutput()
         { if (dataFile != NULL) delete dataFile; }

      protected :

         Poco::HashMap<string, int> getAllLimitedDim(void) const
         {
            Poco::HashMap<string, int> dims;
            const std::set<const CDomain*> sdom = this->getRelFile()->getEnabledDomains();
            const std::set<const CAxis*>  saxis = this->getRelFile()->getEnabledAxis();

            // std::cout << "Nombre de domaines disponibles pour le fichier : " << sdom.size()  << std::endl;
            // std::cout << "Nombre d'axes disponibles pour le fichier : "      << saxis.size() << std::endl;

            std::set<const CDomain*>::const_iterator itt;
            std::set<const CAxis*>::const_iterator it;

            for ( itt = sdom.begin() ; itt != sdom.end(); itt++ )
            {
               string domid = ((*itt)->name.hasValue()) ? (string)(*itt)->name : (*itt)->getId();
               string lonid = string("lon_").append(domid);
               string latid = string("lat_").append(domid);
               dims[lonid]  = (*itt)->ni;
               dims[latid]  = (*itt)->nj;
            }

            for ( it = saxis.begin() ; it != saxis.end(); it++ )
            {
               string axisid = ((*it)->name.hasValue()) ? (string)(*it)->name : (*it)->getId();
               dims[axisid]  = (*it)->size;
            }

            return (dims);
         }

         std::string getFileName(void)
         { return ((getRelFile()->name.hasValue()) ? string(getRelFile()->name): getRelFile()->name.getId()); }

         std::string getTimeStamp(void)
         {
            const int buffer_size = 100;
            time_t rawtime;
            struct tm * timeinfo = NULL;
            char buffer [buffer_size];

            time ( &rawtime );
            timeinfo = localtime ( &rawtime );
            strftime (buffer, buffer_size, "%Y-%b-%d %H:%M:%S %Z", timeinfo);

            return (string(buffer));
         }

         void writeCoords(const string& id, const Array<double, 1>& cdata)
         {
            NcVar *cVar = dataFile->get_var(id.c_str());

            if (!cVar->put(cdata.dataFirst(), cdata.size()))
               throw XMLIOUndefinedValueException("Impossible d'écrire les valeurs de coordonnées "+ id +" !");
         }

      protected : /* virtual */

         virtual void initFile(void)
         {
            string filename = this->getFileName();

            // Création du fichier ou remplacement si celui-ci existe déjà .
            dataFile = new NcFile(filename.c_str(), NcFile::Replace);
            if(!dataFile->is_valid())
               throw XMLIOUndefinedValueException("Impossible d'ouvrir le fichier '"+ filename +"' pour l'écriture des données' !");

            // Ajout de quelques attributs globaux.
            dataFile->add_att("conventions", "CF-1.1");
            dataFile->add_att("file_name"  , getFileName().c_str());
            if (getRelFile()->description.hasValue())
               dataFile->add_att("file_desc" , string(getRelFile()->description).c_str());
            dataFile->add_att("production" , "An IPSL model");
            dataFile->add_att("timeStamp"  , getTimeStamp().c_str()) ;
         }

         virtual void setDims(void)
         {
            bool withTime = true;
            const Poco::HashMap<string, int>& allDim = this->getAllLimitedDim();
            Poco::HashMap<string, int>::ConstIterator it;

            // Ajout des dimensions limitées.
            for (it = allDim.begin() ; it != allDim.end(); it++)
               if (NULL == dataFile->add_dim((*it).first.c_str(), (*it).second))
                  throw XMLIOUndefinedValueException("Impossible d'ajouter la dimension "+ (*it).first +" !");

            // Ajout de la dimension temporelle non limitée.
            if (withTime)
               if (NULL == dataFile->add_dim(TimeName))
                  throw XMLIOUndefinedValueException("Impossible d'ajouter la dimension temporelle !");
         }

         virtual void setCoords(void)
         {
            NcVar *latVar = NULL, // variable de latitude.
                  *lonVar = NULL, // variable de longitude.
                  *othvar = NULL; // toute autre variable.

            Poco::HashMap<string, string> hm;

            const std::set<const CDomain*> sdom = this->getRelFile()->getEnabledDomains();
            const std::set<const CAxis*>  saxis = this->getRelFile()->getEnabledAxis();

            std::set<const CDomain*>::const_iterator itt;
            for ( itt = sdom.begin() ; itt != sdom.end(); itt++ )
            {
               string domid = ((*itt)->name.hasValue()) ? (string)(*itt)->name : (*itt)->getId();
               string lonid = string("lon_").append(domid);
               string latid = string("lat_").append(domid);

               if (!(latVar = dataFile->add_var(latid.c_str(), ncFloat, dataFile->get_dim(latid.c_str()))))
                  throw XMLIOUndefinedValueException("Impossible d'ajouter la variable de latitude !");

               if (!(lonVar = dataFile->add_var(lonid.c_str(), ncFloat, dataFile->get_dim(lonid.c_str()))))
                  throw XMLIOUndefinedValueException("Impossible d'ajouter la variable de longitude !");

               // Attribut de latitude.
               hm["axis"]          = "Y" ;
               hm["standard_name"] = "latitude" ;
               hm["units"]         = "degrees_north";
               hm["long_name"]     = "Latitude" ;
               this->addStringAttributesToVar(latVar, hm);
               hm.clear();

               // Attribut de longitude.
               hm["axis"]           = "X" ;
               hm["standard_name"]  = "longitude" ;
               hm["units"]          = "degrees_east";
               hm["long_name"]      = "Longitude" ;
               this->addStringAttributesToVar(lonVar, hm);
               hm.clear();

               if ((*itt)->yvalue.hasValue())
                  this->writeCoords(lonid, (*itt)->xvalue);
               else throw XMLIOUndefinedValueException("Les coordonnées de longitude (xvalue) ne sont pas définies pour le domaine \""+domid+"\".");

               if ((*itt)->xvalue.hasValue())
                  this->writeCoords(latid, (*itt)->yvalue);
               else throw XMLIOUndefinedValueException("Les coordonnées de latitude (yvalue) ne sont pas définies pour le domaine \""+domid+"\".");
            }

            // Attribut des autres coordonnées.
            std::set<const CAxis*>::const_iterator it;
            for ( it = saxis.begin() ; it != saxis.end(); it++ )
            {
               string axisid = ((*it)->name.hasValue()) ? (string)(*it)->name : (*it)->getId();
               if (!(othvar = dataFile->add_var(axisid.c_str(), ncFloat, dataFile->get_dim(axisid.c_str()))))
                  throw XMLIOUndefinedValueException("Impossible d'ajouter la variable "+ (*it)->getId() +" !");

               hm["axis"] = "Z" ;
               if ((*it)->standard_name.hasValue())hm["standard_name"] = (*it)->standard_name ;
               if ((*it)->long_name.hasValue())    hm["long_name"]     = (*it)->long_name ;
               if ((*it)->unit.hasValue())         hm["units"]         = (*it)->unit;
               this->addStringAttributesToVar(othvar, hm);
               hm.clear();

               if ((*it)->value.hasValue())
                  this->writeCoords(axisid, (*it)->value);
               else throw XMLIOUndefinedValueException("Les coordonnées de l'axe \""+axisid+"\" (value) ne sont pas définies.");
            }
         }

         virtual void setVars(void)
         {
            NcType tvar = ncFloat;
            NcVar *var  = NULL;
            Poco::HashMap<string, string> hm;
            bool lonlat = false, wtime = true;
            const std::vector<CField*>& enabledFields = getRelFile()->getEnabledFields();

            std::vector<CField*>::const_iterator it;

            for ( it = enabledFields.begin() ; it != enabledFields.end(); it++ )
            {
                     CField*  field = (*it);
               const CField* bfield = (*it)->getBaseObject();
               const CDomain*  rdom = field->getGrid()->getRelDomain();
               const CAxis*   raxis = field->getGrid()->getRelAxis();

               string fieldid = (field->name.hasValue()) ? (string)field->name : bfield->getId();
               string domid   = (rdom ->name.hasValue()) ? (string)rdom ->name : rdom  ->getId();
               string axisid  = (raxis->name.hasValue()) ? (string)raxis->name : raxis ->getId();

               string lonid = string("lon_").append(domid); // Nom de la coordonnée longitudinale associée.
               string latid = string("lat_").append(domid); // Nom de la coordonnée latitudinale  associée.

               lonlat = !field->getGrid()->_hasAxis();

               if (field->operation.hasValue())
               { // Si une opération sur le champ est définie ...
                  if (field->operation.getId().compare("once") == 0)
                  { wtime = false; field->freq_op.setValue(NoneDu); }
               }
               else
               { // Si aucune opération sur le champ n'est définie ...
                  FieldOperation ope("inst");
                  field->operation.setValue(ope);
               }

               // Si la fréquence d'opération n'est pas définie, on l'initialise à 0s.
               // Une fréquence à 0s signifie que l'opération se fera à chaque écriture de données depuis la simulation.
               if (!field->freq_op.hasValue())
                  field->freq_op.setValue(NoneDu);
               ((FieldOperation)field->operation).setFreqOp(Duration(field->freq_op));

               tvar = ncFloat;
               if (field->prec.hasValue())
                  if (field->prec == 8) tvar = ncDouble;

               if (wtime)
               {
                  if (lonlat) // 2D spatio + temps
                  {
                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar,
                        dataFile->get_dim(TimeName), dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
                  }
                  else // 3D spatio + temps
                  {
                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar,
                        dataFile->get_dim(TimeName), dataFile->get_dim(axisid.c_str()),
                        dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
                  }
               }
               else
               {
                  if (lonlat) // 2D spatio sans temps
                  {
                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar,
                        dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
                  }
                  else // 3D spatio sans temps
                  {
                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar, dataFile->get_dim(axisid.c_str()),
                        dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
                  }
               }

               if (field->standard_name.hasValue())hm["standard_name"] = field->standard_name ;
               if (field->long_name.hasValue())    hm["long_name"]     = field->long_name ;
               if (field->unit.hasValue())         hm["units"]         = field->unit;
               this->addStringAttributesToVar(var, hm);
               hm.clear();
            }
         }

      private :

         void addStringAttributesToVar(NcVar * const var, const Poco::HashMap<string, string>& attr)
         {
            Poco::HashMap<string, string>::ConstIterator it;

            for ( it = attr.begin() ; it != attr.end(); it++ )
               if (!var->add_att((*it).first.c_str(), (*it).second.c_str()))
                  throw XMLIOUndefinedValueException("Impossible d'ajouter l'attribut' "+ (*it).first +" Ã  la variable "+ var->name() +" !");
         }

         NcFile* dataFile;

   }; //class NetCDF4DataOutput

}// namespace XMLIOSERVER

#endif //__NETCDF4_DATA_OUTPUT__