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

Last change on this file since 128 was 128, checked in by hozdoba, 14 years ago
File size: 15.8 KB
Line 
1#ifndef __NETCDF4_DATA_OUTPUT__
2#define __NETCDF4_DATA_OUTPUT__
3
4#include <netcdfcpp.h>
5
6namespace XMLIOSERVER
7{
8
9   static const char* TimeName = "time";
10
11   class NetCDF4DataOutput : public AbstractDataOutput
12   {
13      public :
14
15         NetCDF4DataOutput(CFile* const _file)
16            : AbstractDataOutput(_file), dataFile(NULL)
17         { /* Ne rien faire de plus */ }
18
19         const NcFile* getDataFile(void) const { return (dataFile); }
20
21         void writeVarData(const string& id, const Array<float, 1>& vdata)
22         {
23            //if (!latVar->put(lats, NLAT))
24            //   return NC_ERR;
25         }
26
27         virtual ~NetCDF4DataOutput()
28         { if (dataFile != NULL) delete dataFile; }
29
30      protected :
31
32         Poco::HashMap<string, int> getAllLimitedDim(void) const
33         {
34            Poco::HashMap<string, int> dims;
35            const std::set<const CDomain*> sdom = this->getRelFile()->getEnabledDomains();
36            const std::set<const CAxis*>  saxis = this->getRelFile()->getEnabledAxis();
37
38            // std::cout << "Nombre de domaines disponibles pour le fichier : " << sdom.size()  << std::endl;
39            // std::cout << "Nombre d'axes disponibles pour le fichier : "      << saxis.size() << std::endl;
40
41            std::set<const CDomain*>::const_iterator itt;
42            std::set<const CAxis*>::const_iterator it;
43
44            for ( itt = sdom.begin() ; itt != sdom.end(); itt++ )
45            {
46               string domid = ((*itt)->name.hasValue()) ? (string)(*itt)->name : (*itt)->getId();
47               string lonid = (sdom.size() == 1)? string("lon"): string("lon_").append(domid);
48               string latid = (sdom.size() == 1)? string("lat"): string("lat_").append(domid);
49               dims[lonid]  = (*itt)->ni;
50               dims[latid]  = (*itt)->nj;
51            }
52
53            for ( it = saxis.begin() ; it != saxis.end(); it++ )
54            {
55               string axisid = ((*it)->name.hasValue()) ? (string)(*it)->name : (*it)->getId();
56               dims[axisid]  = (*it)->size;
57            }
58
59            return (dims);
60         }
61
62         std::string getFileName(void)
63         { return ((getRelFile()->name.hasValue()) ? string(getRelFile()->name): getRelFile()->name.getId()); }
64
65         std::string getTimeStamp(void)
66         {
67            const int buffer_size = 100;
68            time_t rawtime;
69            struct tm * timeinfo = NULL;
70            char buffer [buffer_size];
71
72            time ( &rawtime );
73            timeinfo = localtime ( &rawtime );
74            strftime (buffer, buffer_size, "%Y-%b-%d %H:%M:%S %Z", timeinfo);
75
76            return (string(buffer));
77         }
78
79         void writeCoords(const string& id, const Array<double, 1>& cdata)
80         {
81            NcVar *cVar = dataFile->get_var(id.c_str());
82
83            if (!cVar->put(cdata.dataFirst(), cdata.size()))
84               throw XMLIOUndefinedValueException("Impossible d'écrire les valeurs de coordonnées "+ id +" !");
85         }
86
87      protected : /* virtual */
88
89         virtual void initFile(void)
90         {
91            string filename = this->getFileName();
92
93            // Création du fichier ou remplacement si celui-ci existe déjà.
94            dataFile = new NcFile(filename.c_str(), NcFile::Replace);
95            if(!dataFile->is_valid())
96               throw XMLIOUndefinedValueException("Impossible d'ouvrir le fichier '"+ filename +"' pour l'écriture des données' !");
97
98            // Ajout de quelques attributs globaux.
99            dataFile->add_att("conventions", "CF-1.1");
100            dataFile->add_att("file_name"  , getFileName().c_str());
101            if (getRelFile()->description.hasValue())
102               dataFile->add_att("file_desc" , string(getRelFile()->description).c_str());
103            dataFile->add_att("production" , "An IPSL model");
104            dataFile->add_att("timeStamp"  , getTimeStamp().c_str()) ;
105         }
106
107         virtual void setDims(void)
108         {
109            bool withTime = true;
110            const Poco::HashMap<string, int>& allDim = this->getAllLimitedDim();
111            Poco::HashMap<string, int>::ConstIterator it;
112
113            // Ajout des dimensions limitées.
114            for (it = allDim.begin() ; it != allDim.end(); it++)
115               if (NULL == dataFile->add_dim((*it).first.c_str(), (*it).second))
116                  throw XMLIOUndefinedValueException("Impossible d'ajouter la dimension "+ (*it).first +" !");
117
118            // Ajout de la dimension temporelle non limitée.
119            if (withTime)
120               if (NULL == dataFile->add_dim(TimeName))
121                  throw XMLIOUndefinedValueException("Impossible d'ajouter la dimension temporelle !");
122         }
123
124         virtual void setCoords(void)
125         {
126            NcVar *latVar = NULL, // variable de latitude.
127                  *lonVar = NULL, // variable de longitude.
128                  *othvar = NULL; // toute autre variable.
129
130            Poco::HashMap<string, string> hm;
131            Poco::HashMap<string, float>  hmf;
132
133            const std::set<const CDomain*> sdom = this->getRelFile()->getEnabledDomains();
134            const std::set<const CAxis*>  saxis = this->getRelFile()->getEnabledAxis();
135
136            std::set<const CDomain*>::const_iterator itt;
137            for ( itt = sdom.begin() ; itt != sdom.end(); itt++ )
138            {
139               string domid = ((*itt)->name.hasValue()) ? (string)(*itt)->name : (*itt)->getId();
140               string lonid = (sdom.size() == 1)? string("lon"): string("lon_").append(domid);
141               string latid = (sdom.size() == 1)? string("lat"): string("lat_").append(domid);
142
143               if (!(latVar = dataFile->add_var(latid.c_str(), ncFloat, dataFile->get_dim(latid.c_str()))))
144                  throw XMLIOUndefinedValueException("Impossible d'ajouter la variable de latitude !");
145
146               if (!(lonVar = dataFile->add_var(lonid.c_str(), ncFloat, dataFile->get_dim(lonid.c_str()))))
147                  throw XMLIOUndefinedValueException("Impossible d'ajouter la variable de longitude !");
148
149               // Attribut de latitude.
150               hm["axis"]          = "Y" ;
151               hm["standard_name"] = "latitude" ;
152               hm["units"]         = "degrees_north";
153               hm["long_name"]     = "Latitude" ;
154               hm["nav_model"]     = domid ;
155               if ((*itt)->yvalue.hasValue())
156               {
157                  hmf["valid_min"] = ((Array<double, 1>)(*itt)->yvalue)(0);
158                  hmf["valid_max"] = ((Array<double, 1>)(*itt)->yvalue)(((Array<double, 1>)(*itt)->yvalue).size()-1);
159                  this->addStringAttributesToVar(latVar, hmf);
160                  hmf.clear();
161               }
162               this->addStringAttributesToVar_str(latVar, hm);
163               hm.clear();
164
165               // Attribut de longitude.
166               hm["axis"]          = "X" ;
167               hm["standard_name"] = "longitude" ;
168               hm["units"]         = "degrees_east";
169               hm["long_name"]     = "Longitude" ;
170               hm["nav_model"]     = domid ;
171               if ((*itt)->xvalue.hasValue())
172               {
173                  hmf["valid_min"] = ((Array<double, 1>)(*itt)->xvalue)(0);
174                  hmf["valid_max"] = ((Array<double, 1>)(*itt)->xvalue)(((Array<double, 1>)(*itt)->xvalue).size()-1);
175                  this->addStringAttributesToVar(lonVar, hmf);
176                  hmf.clear();
177               }
178               this->addStringAttributesToVar_str(lonVar, hm);
179               hm.clear();
180
181               if ((*itt)->yvalue.hasValue())
182                  this->writeCoords(lonid, (*itt)->xvalue);
183               else throw XMLIOUndefinedValueException("Les coordonnées de longitude (xvalue) ne sont pas définies pour le domaine \""+domid+"\".");
184
185               if ((*itt)->xvalue.hasValue())
186                  this->writeCoords(latid, (*itt)->yvalue);
187               else throw XMLIOUndefinedValueException("Les coordonnées de latitude (yvalue) ne sont pas définies pour le domaine \""+domid+"\".");
188            }
189
190            // Attribut des autres coordonnées.
191            std::set<const CAxis*>::const_iterator it;
192            for ( it = saxis.begin() ; it != saxis.end(); it++ )
193            {
194               string axisid = ((*it)->name.hasValue()) ? (string)(*it)->name : (*it)->getId();
195               if (!(othvar = dataFile->add_var(axisid.c_str(), ncFloat, dataFile->get_dim(axisid.c_str()))))
196                  throw XMLIOUndefinedValueException("Impossible d'ajouter la variable "+ (*it)->getId() +" !");
197
198               hm["axis"]     = "Z" ;
199               hm["title"]    = axisid ;
200               hm["positive"] = "unknown" ;
201
202               if ((*it)->standard_name.hasValue())hm["standard_name"] = (*it)->standard_name ;
203               if ((*it)->long_name.hasValue())    hm["long_name"]     = (*it)->long_name ;
204               if ((*it)->unit.hasValue())         hm["units"]         = (*it)->unit;
205               if ((*it)->value.hasValue())
206               {
207                  hmf["valid_min"] = ((Array<double, 1>)(*it)->value)(0);
208                  hmf["valid_max"] = ((Array<double, 1>)(*it)->value)(((Array<double, 1>)(*it)->value).size()-1);
209                  this->addStringAttributesToVar(othvar, hmf);
210                  hmf.clear();
211               }
212               this->addStringAttributesToVar_str(othvar, hm);
213               hm.clear();
214
215               if ((*it)->value.hasValue())
216                  this->writeCoords(axisid, (*it)->value);
217               else throw XMLIOUndefinedValueException("Les coordonnées de l'axe \""+axisid+"\" (value) ne sont pas définies.");
218            }
219         }
220
221         virtual void setVars(void)
222         {
223            NcType tvar = ncFloat;
224            NcVar *var  = NULL;
225            Poco::HashMap<string, string> hm;
226            Poco::HashMap<string, float>  hmf;
227            bool lonlat = false, wtime = true;
228            const std::vector<CField*>& enabledFields = getRelFile()->getEnabledFields();
229            const std::set<const CDomain*> sdom = this->getRelFile()->getEnabledDomains();
230
231            std::vector<CField*>::const_iterator it;
232
233            for ( it = enabledFields.begin() ; it != enabledFields.end(); it++ )
234            {
235                     CField  *  field = (*it);
236               const CField  * bfield = (*it)->getBaseObject();
237               const CDomain *   rdom = field->getGrid()->getRelDomain();
238               const CAxis   *  raxis = field->getGrid()->getRelAxis();
239               FieldOperation*   fope = field->operation.getValue();
240               Duration      *  ffope = field->freq_op.getValue();
241
242               string fieldid = (field->name.hasValue()) ? (string)field->name : bfield->getId();
243               string domid   = (rdom ->name.hasValue()) ? (string)rdom ->name : rdom  ->getId();
244               string axisid  = (raxis->name.hasValue()) ? (string)raxis->name : raxis ->getId();
245
246               string lonid = (sdom.size() == 1)? string("lon"): string("lon_").append(domid); // Nom de la coordonnée longitudinale associée.
247               string latid = (sdom.size() == 1)? string("lat"): string("lat_").append(domid); // Nom de la coordonnée latitudinale  associée.
248
249               lonlat = !field->getGrid()->_hasAxis();
250
251
252               ////////////////////////////////////// A DEPLACER
253               if (fope != NULL)
254               { // Si une opération sur le champ est définie ...
255                  if (fope->getId().compare("once") == 0)
256                  { wtime = false; field->freq_op.setValue(NoneDu); }
257               }
258               else
259               { // Si aucune opération sur le champ n'est définie ...
260                  FieldOperation ope("inst");
261                  field->operation.setValue(ope);
262               }
263
264               // Si la fréquence d'opération n'est pas définie, on l'initialise à 0s.
265               // Une fréquence à 0s signifie que l'opération se fera à chaque écriture de données depuis la simulation.
266               if (ffope == NULL)
267                  field->freq_op.setValue(NoneDu);
268
269               field->operation.getValue()->setFreqOp(field->freq_op);
270               /////////////////////////////////////////////
271
272               tvar = ncFloat;
273               if (field->prec.hasValue())
274                  if (field->prec == 8) tvar = ncDouble;
275
276               if (wtime)
277               {
278                  if (lonlat) // 2D spatio + temps
279                  {
280                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar,
281                        dataFile->get_dim(TimeName), dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
282                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
283                  }
284                  else // 3D spatio + temps
285                  {
286                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar,
287                        dataFile->get_dim(TimeName), dataFile->get_dim(axisid.c_str()),
288                        dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
289                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
290                  }
291               }
292               else
293               {
294                  if (lonlat) // 2D spatio sans temps
295                  {
296                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar,
297                        dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
298                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
299                  }
300                  else // 3D spatio sans temps
301                  {
302                     if (!(var = dataFile->add_var(fieldid.c_str(), tvar, dataFile->get_dim(axisid.c_str()),
303                        dataFile->get_dim(latid.c_str()), dataFile->get_dim(lonid.c_str()))))
304                        throw XMLIOUndefinedValueException("Impossible d'ajouter le champ "+ field->getId() +" !");
305                  }
306               }
307
308               if (field->standard_name.hasValue())hm["standard_name"]     = field->standard_name ;
309               if (field->long_name.hasValue())    hm["long_name"]         = field->long_name ;
310               if (field->unit.hasValue())         hm["units"]             = field->unit;
311
312               if (field->operation.hasValue())
313               {
314                  hm["online_operation"]   = field->operation.getValue()->getId();
315
316                  if (!field->operation.getValue()->getFreqOp().isNone())
317                     hm["interval_operation"] = field->operation.getValue()->getFreqOp().toString();
318                  else
319                     hm["interval_operation"] = ((Duration)this->getRelFile()->output_freq).toString();
320               }
321               hm["interval_write"] = ((Duration)this->getRelFile()->output_freq).toString();
322
323               this->addStringAttributesToVar_str(var, hm);
324               hm.clear();
325            }
326         }
327
328      private :
329
330         template <class U>
331            void addStringAttributesToVar(NcVar * const var, const Poco::HashMap<string, U>& attr)
332         {
333            typename Poco::HashMap<string, U>::ConstIterator it;
334
335            for ( it = attr.begin() ; it != attr.end(); it++ )
336               if (!var->add_att((*it).first.c_str(), (*it).second))
337                  throw XMLIOUndefinedValueException("Impossible d'ajouter l'attribut' "+ (*it).first +" à la variable "+ var->name() +" !");
338         }
339
340         void addStringAttributesToVar_str(NcVar * const var, const Poco::HashMap<string, string>& attr)
341         {
342            Poco::HashMap<string, string>::ConstIterator it;
343
344            for ( it = attr.begin() ; it != attr.end(); it++ )
345               if (!var->add_att((*it).first.c_str(), (*it).second.c_str()))
346                  throw XMLIOUndefinedValueException("Impossible d'ajouter l'attribut' "+ (*it).first +" à la variable "+ var->name() +" !");
347         }
348
349         NcFile* dataFile;
350
351   }; //class NetCDF4DataOutput
352
353}// namespace XMLIOSERVER
354
355#endif //__NETCDF4_DATA_OUTPUT__
Note: See TracBrowser for help on using the repository browser.