source: XIOS/dev/branch_yushan/src/node/domain.cpp @ 1109

Last change on this file since 1109 was 1109, checked in by yushan, 7 years ago

test_omp OK

  • Property copyright set to
    Software name : XIOS (Xml I/O Server)
    http://forge.ipsl.jussieu.fr/ioserver
    Creation date : January 2009
    Licence : CeCCIL version2
    see license file in root directory : Licence_CeCILL_V2-en.txt
    or http://www.cecill.info/licences/Licence_CeCILL_V2-en.html
    Holder : CEA/LSCE (Laboratoire des Sciences du CLimat et de l'Environnement)
    CNRS/IPSL (Institut Pierre Simon Laplace)
    Project Manager : Yann Meurdesoif
    yann.meurdesoif@cea.fr
  • Property svn:executable set to *
File size: 78.1 KB
Line 
1#include "domain.hpp"
2
3#include "attribute_template.hpp"
4#include "object_template.hpp"
5#include "group_template.hpp"
6
7#include "xios_spl.hpp"
8#include "event_client.hpp"
9#include "event_server.hpp"
10#include "buffer_in.hpp"
11#include "message.hpp"
12#include "type.hpp"
13#include "context.hpp"
14#include "context_client.hpp"
15#include "context_server.hpp"
16#include "array_new.hpp"
17#include "distribution_client.hpp"
18#include "server_distribution_description.hpp"
19#include "client_server_mapping_distributed.hpp"
20#include "zoom_domain.hpp"
21#include "interpolate_domain.hpp"
22#include "generate_rectilinear_domain.hpp"
23
24#include <algorithm>
25
26namespace xios {
27
28   /// ////////////////////// Définitions ////////////////////// ///
29
30   CDomain::CDomain(void)
31      : CObjectTemplate<CDomain>(), CDomainAttributes()
32      , isChecked(false), relFiles(), isClientChecked(false), nbConnectedClients_(), indSrv_(), connectedServerRank_()
33      , hasBounds(false), hasArea(false), isDistributed_(false), nGlobDomain_(), isCompressible_(false), isUnstructed_(false)
34      , isClientAfterTransformationChecked(false), hasLonLat(false)
35      , lonvalue_client(), latvalue_client(), bounds_lon_client(), bounds_lat_client()
36      , isRedistributed_(false), hasPole(false)
37   {
38   }
39
40   CDomain::CDomain(const StdString & id)
41      : CObjectTemplate<CDomain>(id), CDomainAttributes()
42      , isChecked(false), relFiles(), isClientChecked(false), nbConnectedClients_(), indSrv_(), connectedServerRank_()
43      , hasBounds(false), hasArea(false), isDistributed_(false), nGlobDomain_(), isCompressible_(false), isUnstructed_(false)
44      , isClientAfterTransformationChecked(false), hasLonLat(false)
45      , lonvalue_client(), latvalue_client(), bounds_lon_client(), bounds_lat_client()
46      , isRedistributed_(false), hasPole(false)
47   {
48         }
49
50   CDomain::~CDomain(void)
51   {
52   }
53
54   ///---------------------------------------------------------------
55
56   void CDomain::assignMesh(const StdString meshName, const int nvertex)
57   {
58     mesh = CMesh::getMesh(meshName, nvertex);
59   }
60
61   CDomain* CDomain::createDomain()
62   {
63     CDomain* domain = CDomainGroup::get("domain_definition")->createChild();
64     return domain;
65   }
66
67   //std::map<StdString, ETranformationType> CDomain::transformationMapList_ = std::map<StdString, ETranformationType>();
68   //bool CDomain::_dummyTransformationMapList = CDomain::initializeTransformationMap(CDomain::transformationMapList_);
69
70   std::map<StdString, ETranformationType> *CDomain::transformationMapList_ptr = 0;
71
72   bool CDomain::initializeTransformationMap(std::map<StdString, ETranformationType>& m)
73   {
74     m["zoom_domain"] = TRANS_ZOOM_DOMAIN;
75     m["interpolate_domain"] = TRANS_INTERPOLATE_DOMAIN;
76     m["generate_rectilinear_domain"] = TRANS_GENERATE_RECTILINEAR_DOMAIN;
77     m["compute_connectivity_domain"] = TRANS_COMPUTE_CONNECTIVITY_DOMAIN;
78     m["expand_domain"] = TRANS_EXPAND_DOMAIN;
79   }
80
81   bool CDomain::initializeTransformationMap()
82   {
83     CDomain::transformationMapList_ptr = new std::map<StdString, ETranformationType>();
84     (*CDomain::transformationMapList_ptr)["zoom_domain"] = TRANS_ZOOM_DOMAIN;
85     (*CDomain::transformationMapList_ptr)["interpolate_domain"] = TRANS_INTERPOLATE_DOMAIN;
86     (*CDomain::transformationMapList_ptr)["generate_rectilinear_domain"] = TRANS_GENERATE_RECTILINEAR_DOMAIN;
87     (*CDomain::transformationMapList_ptr)["compute_connectivity_domain"] = TRANS_COMPUTE_CONNECTIVITY_DOMAIN;
88     (*CDomain::transformationMapList_ptr)["expand_domain"] = TRANS_EXPAND_DOMAIN;
89   }
90
91
92   const std::set<StdString> & CDomain::getRelFiles(void) const
93   {
94      return (this->relFiles);
95   }
96
97
98   const std::vector<int>& CDomain::getIndexesToWrite(void) const
99   {
100     return indexesToWrite;
101   }
102
103   /*!
104     Returns the number of indexes written by each server.
105     \return the number of indexes written by each server
106   */
107   int CDomain::getNumberWrittenIndexes() const
108   {
109     return numberWrittenIndexes_;
110   }
111
112   /*!
113     Returns the total number of indexes written by the servers.
114     \return the total number of indexes written by the servers
115   */
116   int CDomain::getTotalNumberWrittenIndexes() const
117   {
118     return totalNumberWrittenIndexes_;
119   }
120
121   /*!
122     Returns the offset of indexes written by each server.
123     \return the offset of indexes written by each server
124   */
125   int CDomain::getOffsetWrittenIndexes() const
126   {
127     return offsetWrittenIndexes_;
128   }
129
130   //----------------------------------------------------------------
131
132   /*!
133    * Compute the minimum buffer size required to send the attributes to the server(s).
134    *
135    * \return A map associating the server rank with its minimum buffer size.
136    */
137   std::map<int, StdSize> CDomain::getAttributesBufferSize()
138   {
139     CContextClient* client = CContext::getCurrent()->client;
140
141     std::map<int, StdSize> attributesSizes = getMinimumBufferSizeForAttributes();
142
143     if (client->isServerLeader())
144     {
145       // size estimation for sendServerAttribut
146       size_t size = 11 * sizeof(size_t);
147
148       const std::list<int>& ranks = client->getRanksServerLeader();
149       for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
150       {
151         if (size > attributesSizes[*itRank])
152           attributesSizes[*itRank] = size;
153       }
154     }
155
156     std::map<int, std::vector<size_t> >::const_iterator itIndexEnd = indSrv_.end();
157     std::map<int, std::vector<int> >::const_iterator itWrittenIndexEnd = indWrittenSrv_.end();
158     for (size_t k = 0; k < connectedServerRank_.size(); ++k)
159     {
160       int rank = connectedServerRank_[k];
161       std::map<int, std::vector<size_t> >::const_iterator it = indSrv_.find(rank);
162       size_t idxCount = (it != itIndexEnd) ? it->second.size() : 0;
163
164       // size estimation for sendIndex (and sendArea which is always smaller or equal)
165       size_t sizeIndexEvent = 2 * sizeof(size_t) + 2 * CArray<int,1>::size(idxCount);
166       if (isCompressible_)
167       {
168         std::map<int, std::vector<int> >::const_iterator itWritten = indWrittenSrv_.find(rank);
169         size_t writtenIdxCount = (itWritten != itWrittenIndexEnd) ? itWritten->second.size() : 0;
170         sizeIndexEvent += CArray<int,1>::size(writtenIdxCount);
171       }
172
173       // size estimation for sendLonLat
174       size_t sizeLonLatEvent = CArray<double,1>::size(idxCount);
175       if (hasBounds)
176         sizeLonLatEvent += CArray<double,2>::size(nvertex * idxCount);
177
178       size_t size = CEventClient::headerSize + getId().size() + sizeof(size_t) + std::max(sizeIndexEvent, sizeLonLatEvent);
179       if (size > attributesSizes[rank])
180         attributesSizes[rank] = size;
181     }
182
183     return attributesSizes;
184   }
185
186   //----------------------------------------------------------------
187
188   bool CDomain::isEmpty(void) const
189   {
190      return ((this->zoom_ni_srv == 0) ||
191              (this->zoom_nj_srv == 0));
192   }
193
194   //----------------------------------------------------------------
195
196   bool CDomain::IsWritten(const StdString & filename) const
197   {
198      return (this->relFiles.find(filename) != this->relFiles.end());
199   }
200
201   bool CDomain::isWrittenCompressed(const StdString& filename) const
202   {
203      return (this->relFilesCompressed.find(filename) != this->relFilesCompressed.end());
204   }
205
206   //----------------------------------------------------------------
207
208   bool CDomain::isDistributed(void) const
209   {
210      return isDistributed_;
211   }
212
213   //----------------------------------------------------------------
214
215   /*!
216    * Test whether the data defined on the domain can be outputted in a compressed way.
217    *
218    * \return true if and only if a mask was defined for this domain
219    */
220   bool CDomain::isCompressible(void) const
221   {
222      return isCompressible_;
223   }
224
225   void CDomain::addRelFile(const StdString & filename)
226   {
227      this->relFiles.insert(filename);
228   }
229
230   void CDomain::addRelFileCompressed(const StdString& filename)
231   {
232      this->relFilesCompressed.insert(filename);
233   }
234
235   StdString CDomain::GetName(void)   { return (StdString("domain")); }
236   StdString CDomain::GetDefName(void){ return (CDomain::GetName()); }
237   ENodeType CDomain::GetType(void)   { return (eDomain); }
238
239   //----------------------------------------------------------------
240
241   /*!
242     Redistribute RECTILINEAR domain with a number of local domains.
243   All attributes ni,nj,ibegin,jbegin (if defined) will be rewritten
244   The optional attributes lonvalue, latvalue will be added. Because this function only serves (for now)
245   for interpolation from unstructured domain to rectilinear one, range of latvalue is 0-360 and lonvalue is -90 - +90
246    \param [in] nbLocalDomain number of local domain on the domain destination
247   */
248   void CDomain::redistribute(int nbLocalDomain)
249   {
250     if (this->isRedistributed_) return;
251
252     this->isRedistributed_ = true;
253     CContext* context = CContext::getCurrent();
254     CContextClient* client = context->client;
255     int rankClient = client->clientRank;
256     int rankOnDomain = rankClient%nbLocalDomain;
257
258     if (ni_glo.isEmpty() || ni_glo <= 0 )
259     {
260        ERROR("CDomain::redistribute(int nbLocalDomain)",
261           << "[ Id = " << this->getId() << " ] "
262           << "The global domain is badly defined,"
263           << " check the \'ni_glo\'  value !")
264     }
265
266     if (nj_glo.isEmpty() || nj_glo <= 0 )
267     {
268        ERROR("CDomain::redistribute(int nbLocalDomain)",
269           << "[ Id = " << this->getId() << " ] "
270           << "The global domain is badly defined,"
271           << " check the \'nj_glo\'  value !")
272     }
273
274     if ((type_attr::rectilinear == type)  || (type_attr::curvilinear == type))
275     {
276        int globalDomainSize = ni_glo * nj_glo;
277        if (globalDomainSize <= nbLocalDomain)
278        {
279          for (int idx = 0; idx < nbLocalDomain; ++idx)
280          {
281            if (rankOnDomain < globalDomainSize)
282            {
283              int iIdx = rankOnDomain % ni_glo;
284              int jIdx = rankOnDomain / ni_glo;
285              ibegin.setValue(iIdx); jbegin.setValue(jIdx);
286              ni.setValue(1); nj.setValue(1);
287            }
288            else
289            {
290              ibegin.setValue(0); jbegin.setValue(0);
291              ni.setValue(0); nj.setValue(0);
292            }
293          }
294        }
295        else
296        {
297          float njGlo = nj_glo.getValue();
298          float niGlo = ni_glo.getValue();
299          int nbProcOnX, nbProcOnY, range;
300
301          // Compute (approximately) number of segment on x and y axis
302          float yOverXRatio = njGlo/niGlo;
303
304          nbProcOnX = std::ceil(std::sqrt(nbLocalDomain/yOverXRatio));
305          nbProcOnY = std::ceil(((float)nbLocalDomain)/nbProcOnX);
306
307          // Simple distribution: Sweep from top to bottom, left to right
308          // Calculate local begin on x
309          std::vector<int> ibeginVec(nbProcOnX,0), jbeginVec(nbProcOnY,0);
310          std::vector<int> niVec(nbProcOnX), njVec(nbProcOnY);
311          for (int i = 1; i < nbProcOnX; ++i)
312          {
313            range = ni_glo / nbProcOnX;
314            if (i < (ni_glo%nbProcOnX)) ++range;
315            niVec[i-1] = range;
316            ibeginVec[i] = ibeginVec[i-1] + niVec[i-1];
317          }
318          niVec[nbProcOnX-1] = ni_glo - ibeginVec[nbProcOnX-1];
319
320          // Calculate local begin on y
321          for (int j = 1; j < nbProcOnY; ++j)
322          {
323            range = nj_glo / nbProcOnY;
324            if (j < (nj_glo%nbProcOnY)) ++range;
325            njVec[j-1] = range;
326            jbeginVec[j] = jbeginVec[j-1] + njVec[j-1];
327          }
328          njVec[nbProcOnY-1] = nj_glo - jbeginVec[nbProcOnY-1];
329
330          // Now assign value to ni, ibegin, nj, jbegin
331          int iIdx = rankOnDomain % nbProcOnX;
332          int jIdx = rankOnDomain / nbProcOnX;
333
334          if (rankOnDomain != (nbLocalDomain-1))
335          {
336            ibegin.setValue(ibeginVec[iIdx]);
337            jbegin.setValue(jbeginVec[jIdx]);
338            nj.setValue(njVec[jIdx]);
339            ni.setValue(niVec[iIdx]);
340          }
341          else // just merge all the remaining rectangle into the last one
342          {
343            ibegin.setValue(ibeginVec[iIdx]);
344            jbegin.setValue(jbeginVec[jIdx]);
345            nj.setValue(njVec[jIdx]);
346            ni.setValue(ni_glo - ibeginVec[iIdx]);
347          }
348        }
349
350        // Now fill other attributes
351        if (type_attr::rectilinear == type) fillInRectilinearLonLat();
352     }
353     else  // unstructured domain
354     {
355       if (this->i_index.isEmpty())
356       {
357          int globalDomainSize = ni_glo * nj_glo;
358          if (globalDomainSize <= nbLocalDomain)
359          {
360            for (int idx = 0; idx < nbLocalDomain; ++idx)
361            {
362              if (rankOnDomain < globalDomainSize)
363              {
364                int iIdx = rankOnDomain % ni_glo;
365                int jIdx = rankOnDomain / ni_glo;
366                ibegin.setValue(iIdx); jbegin.setValue(jIdx);
367                ni.setValue(1); nj.setValue(1);
368              }
369              else
370              {
371                ibegin.setValue(0); jbegin.setValue(0);
372                ni.setValue(0); nj.setValue(0);
373              }
374            }
375          }
376          else
377          {
378            float njGlo = nj_glo.getValue();
379            float niGlo = ni_glo.getValue();
380            std::vector<int> ibeginVec(nbLocalDomain,0);
381            std::vector<int> niVec(nbLocalDomain);
382            for (int i = 1; i < nbLocalDomain; ++i)
383            {
384              int range = ni_glo / nbLocalDomain;
385              if (i < (ni_glo%nbLocalDomain)) ++range;
386              niVec[i-1] = range;
387              ibeginVec[i] = ibeginVec[i-1] + niVec[i-1];
388            }
389            niVec[nbLocalDomain-1] = ni_glo - ibeginVec[nbLocalDomain-1];
390
391            int iIdx = rankOnDomain % nbLocalDomain;
392            ibegin.setValue(ibeginVec[iIdx]);
393            jbegin.setValue(0);
394            ni.setValue(niVec[iIdx]);
395            nj.setValue(1);
396          }
397        }
398        else
399        {
400          ibegin.setValue(this->i_index(0));
401          jbegin.setValue(0);
402          ni.setValue(this->i_index.numElements());
403          nj.setValue(1);
404        }
405     }
406
407     checkDomain();
408   }
409
410   /*!
411     Fill in the values for lonvalue_1d and latvalue_1d of rectilinear domain
412     Range of longitude value from 0 - 360
413     Range of latitude value from -90 - +90
414   */
415   void CDomain::fillInRectilinearLonLat()
416   {
417     if (!lonvalue_rectilinear_read_from_file.isEmpty())
418     {
419       lonvalue_1d.resize(ni);
420       for (int idx = 0; idx < ni; ++idx)
421         lonvalue_1d(idx) = lonvalue_rectilinear_read_from_file(idx+ibegin);
422       lon_start.setValue(lonvalue_rectilinear_read_from_file(0));
423       lon_end.setValue(lonvalue_rectilinear_read_from_file(ni_glo-1));
424     }
425     else
426     {
427       if (!lonvalue_2d.isEmpty()) lonvalue_2d.free();
428       lonvalue_1d.resize(ni);
429       double lonRange = lon_end - lon_start;
430       double lonStep = (1 == ni_glo.getValue()) ? lonRange : lonRange/double(ni_glo.getValue()-1);
431
432        // Assign lon value
433       for (int i = 0; i < ni; ++i)
434       {
435         if (0 == (ibegin + i))
436         {
437           lonvalue_1d(i) = lon_start;
438         }
439         else if (ni_glo == (ibegin + i + 1))
440         {
441           lonvalue_1d(i) = lon_end;
442         }
443         else
444         {
445           lonvalue_1d(i) = (ibegin + i) * lonStep  + lon_start;
446         }
447       }
448     }
449
450
451     if (!latvalue_rectilinear_read_from_file.isEmpty())
452     {
453       latvalue_1d.resize(nj);
454       for (int idx = 0; idx < nj; ++idx)
455         latvalue_1d(idx) = latvalue_rectilinear_read_from_file(idx+jbegin);
456       lat_start.setValue(latvalue_rectilinear_read_from_file(0));
457       lat_end.setValue(latvalue_rectilinear_read_from_file(nj_glo-1));
458     }
459     else
460     {
461       if (!latvalue_2d.isEmpty()) latvalue_1d.free();
462       latvalue_1d.resize(nj);
463
464       double latRange = lat_end - lat_start;
465       double latStep = (1 == nj_glo.getValue()) ? latRange : latRange/double(nj_glo.getValue()-1);
466
467       for (int j = 0; j < nj; ++j)
468       {
469         if (0 == (jbegin + j))
470         {
471            latvalue_1d(j) = lat_start;
472         }
473         else if (nj_glo == (jbegin + j + 1))
474         {
475            latvalue_1d(j) = lat_end;
476         }
477         else
478         {
479           latvalue_1d(j) =  (jbegin + j) * latStep + lat_start;
480         }
481       }
482     }
483   }
484
485
486
487   void CDomain::AllgatherRectilinearLonLat(CArray<double,1>& lon, CArray<double,1>& lat, CArray<double,1>& lon_g, CArray<double,1>& lat_g)
488   {
489          CContext* context = CContext::getCurrent();
490    CContextClient* client = context->client;
491          lon_g.resize(ni_glo) ;
492          lat_g.resize(nj_glo) ;
493
494
495          int* ibegin_g = new int[client->clientSize] ;
496          int* jbegin_g = new int[client->clientSize] ;
497          int* ni_g = new int[client->clientSize] ;
498          int* nj_g = new int[client->clientSize] ;
499          int v ;
500          v=ibegin ;
501          MPI_Allgather(&v,1,MPI_INT,ibegin_g,1,MPI_INT,client->intraComm) ;
502          v=jbegin ;
503          MPI_Allgather(&v,1,MPI_INT,jbegin_g,1,MPI_INT,client->intraComm) ;
504          v=ni ;
505          MPI_Allgather(&v,1,MPI_INT,ni_g,1,MPI_INT,client->intraComm) ;
506          v=nj ;
507          MPI_Allgather(&v,1,MPI_INT,nj_g,1,MPI_INT,client->intraComm) ;
508
509          MPI_Allgatherv(lon.dataFirst(),ni,MPI_DOUBLE,lon_g.dataFirst(),ni_g, ibegin_g,MPI_DOUBLE,client->intraComm) ;
510          MPI_Allgatherv(lat.dataFirst(),nj,MPI_DOUBLE,lat_g.dataFirst(),nj_g, jbegin_g,MPI_DOUBLE,client->intraComm) ;
511
512      delete[] ibegin_g ;
513      delete[] jbegin_g ;
514      delete[] ni_g ;
515      delete[] nj_g ;
516   }
517
518   void CDomain::fillInRectilinearBoundLonLat(CArray<double,1>& lon, CArray<double,1>& lat,
519                                              CArray<double,2>& boundsLon, CArray<double,2>& boundsLat)
520   {
521     int i,j,k;
522
523     const int nvertexValue = 4;
524     boundsLon.resize(nvertexValue,ni*nj);
525
526     if (ni_glo>1)
527     {
528       double lonStepStart = lon(1)-lon(0);
529       bounds_lon_start=lon(0) - lonStepStart/2;
530       double lonStepEnd = lon(ni_glo-1)-lon(ni_glo-2);
531       bounds_lon_end=lon(ni_glo-1) + lonStepEnd/2;
532       double errorBoundsLon = std::abs(360-std::abs(bounds_lon_end-bounds_lon_start));
533
534       // if errorBoundsLon is reasonably small (0.1 x cell size) consider it as closed in longitude
535       if (errorBoundsLon < std::abs(lonStepStart)*1e-1 || errorBoundsLon < std::abs(lonStepEnd)*1e-1 )
536       {
537         bounds_lon_start= (lon(0) + lon(ni_glo-1)-360)/2 ;
538         bounds_lon_end= (lon(0) +360 + lon(ni_glo-1))/2 ;
539       }
540     }
541     else
542     {
543       if (bounds_lon_start.isEmpty()) bounds_lon_start=-180. ;
544       if (bounds_lon_end.isEmpty()) bounds_lon_end=180.-1e-8 ;
545     }
546
547     for(j=0;j<nj;++j)
548       for(i=0;i<ni;++i)
549       {
550         k=j*ni+i;
551         boundsLon(0,k) = boundsLon(1,k) = (0 == (ibegin + i)) ? bounds_lon_start
552                                                               : (lon(ibegin + i)+lon(ibegin + i-1))/2;
553         boundsLon(2,k) = boundsLon(3,k) = ((ibegin + i + 1) == ni_glo) ? bounds_lon_end
554                                                                        : (lon(ibegin + i + 1)+lon(ibegin + i))/2;
555       }
556
557
558    boundsLat.resize(nvertexValue,nj*ni);
559    bool isNorthPole=false ;
560    bool isSouthPole=false ;
561    if (std::abs(90 - std::abs(lat(0))) < NumTraits<double>::epsilon()) isNorthPole = true;
562    if (std::abs(90 - std::abs(lat(nj_glo-1))) < NumTraits<double>::epsilon()) isSouthPole = true;
563
564    // lat boundaries beyond pole the assimilate it to pole
565    // lat boundarie is relativelly close to pole (0.1 x cell size) assimilate it to pole
566    if (nj_glo>1)
567    {
568      double latStepStart = lat(1)-lat(0);
569      if (isNorthPole) bounds_lat_start=lat(0);
570      else
571      {
572        bounds_lat_start=lat(0)-latStepStart/2;
573        if (bounds_lat_start >= 90 ) bounds_lat_start=90 ;
574        else if (bounds_lat_start <= -90 ) bounds_lat_start=-90 ;
575        else if (bounds_lat_start <= 90 && bounds_lat_start >= lat(0))
576        {
577          if ( std::abs(90-bounds_lat_start) <= 0.1*std::abs(latStepStart)) bounds_lat_start=90 ;
578        }
579        else if (bounds_lat_start >= -90 && bounds_lat_start <= lat(0))
580        {
581          if ( std::abs(-90 - bounds_lat_start) <= 0.1*std::abs(latStepStart)) bounds_lat_start=-90 ;
582        }
583      }
584
585      double latStepEnd = lat(nj_glo-1)-lat(nj_glo-2);
586      if (isSouthPole) bounds_lat_end=lat(nj_glo-1);
587      else
588      {
589        bounds_lat_end=lat(nj_glo-1)+latStepEnd/2;
590
591        if (bounds_lat_end >= 90 ) bounds_lat_end=90 ;
592        else if (bounds_lat_end <= -90 ) bounds_lat_end=-90 ;
593        else if (bounds_lat_end <= 90 && bounds_lat_end >= lat(nj_glo-1))
594        {
595          if ( std::abs(90-bounds_lat_end) <= 0.1*std::abs(latStepEnd)) bounds_lat_end=90 ;
596        }
597        else if (bounds_lat_end >= -90 && bounds_lat_end <= lat(nj_glo-1))
598        {
599          if ( std::abs(-90 - bounds_lat_end) <= 0.1*std::abs(latStepEnd)) bounds_lat_end=-90 ;
600        }
601      }
602    }
603    else
604    {
605      if (bounds_lat_start.isEmpty()) bounds_lon_start=-90. ;
606      if (bounds_lat_end.isEmpty()) bounds_lat_end=90 ;
607    }
608
609    for(j=0;j<nj;++j)
610      for(i=0;i<ni;++i)
611      {
612        k=j*ni+i;
613        boundsLat(1,k) = boundsLat(2,k) = (0 == (jbegin + j)) ? bounds_lat_start
614                                                              : (lat(jbegin + j)+lat(jbegin + j-1))/2;
615        boundsLat(0,k) = boundsLat(3,k) = ((jbegin + j +1) == nj_glo) ? bounds_lat_end
616                                                                      : (lat(jbegin + j + 1)+lat(jbegin + j))/2;
617      }
618   }
619
620   void CDomain::checkDomain(void)
621   {
622     if (type.isEmpty())
623     {
624       ERROR("CDomain::checkDomain(void)",
625             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
626             << "The domain type is mandatory, "
627             << "please define the 'type' attribute.")
628     }
629
630     if (type == type_attr::gaussian) 
631     {
632           hasPole=true ;
633             type.setValue(type_attr::unstructured) ;
634           }
635           else if (type == type_attr::rectilinear) hasPole=true ;
636         
637     if (type == type_attr::unstructured)
638     {
639        if (ni_glo.isEmpty())
640        {
641          ERROR("CDomain::checkDomain(void)",
642                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
643                << "The global domain is badly defined, "
644                << "the mandatory 'ni_glo' attribute is missing.")
645        }
646        else if (ni_glo <= 0)
647        {
648          ERROR("CDomain::checkDomain(void)",
649                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
650                << "The global domain is badly defined, "
651                << "'ni_glo' attribute should be strictly positive so 'ni_glo = " << ni_glo.getValue() << "' is invalid.")
652        }
653        isUnstructed_ = true;
654        nj_glo = 1;
655        nj = 1;
656        jbegin = 0;
657        if (!i_index.isEmpty()) ni = i_index.numElements();
658        j_index.resize(ni);
659        for(int i=0;i<ni;++i) j_index(i)=0;
660
661        if (!area.isEmpty())
662          area.transposeSelf(1, 0);
663     }
664
665     if (ni_glo.isEmpty())
666     {
667       ERROR("CDomain::checkDomain(void)",
668             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
669             << "The global domain is badly defined, "
670             << "the mandatory 'ni_glo' attribute is missing.")
671     }
672     else if (ni_glo <= 0)
673     {
674       ERROR("CDomain::checkDomain(void)",
675             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
676             << "The global domain is badly defined, "
677             << "'ni_glo' attribute should be strictly positive so 'ni_glo = " << ni_glo.getValue() << "' is invalid.")
678     }
679
680     if (nj_glo.isEmpty())
681     {
682       ERROR("CDomain::checkDomain(void)",
683             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
684             << "The global domain is badly defined, "
685             << "the mandatory 'nj_glo' attribute is missing.")
686     }
687     else if (nj_glo <= 0)
688     {
689       ERROR("CDomain::checkDomain(void)",
690             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
691             << "The global domain is badly defined, "
692             << "'nj_glo' attribute should be strictly positive so 'nj_glo = " << nj_glo.getValue() << "' is invalid.")
693     }
694
695     checkLocalIDomain();
696     checkLocalJDomain();
697
698     if (i_index.isEmpty())
699     {
700       i_index.resize(ni*nj);
701       for (int j = 0; j < nj; ++j)
702         for (int i = 0; i < ni; ++i) i_index(i+j*ni) = i+ibegin;
703     }
704
705     if (j_index.isEmpty())
706     {
707       j_index.resize(ni*nj);
708       for (int j = 0; j < nj; ++j)
709         for (int i = 0; i < ni; ++i) j_index(i+j*ni) = j+jbegin;
710     }
711     computeNGlobDomain();
712     checkZoom();
713
714     isDistributed_ = !((!ni.isEmpty() && (ni == ni_glo) && !nj.isEmpty() && (nj == nj_glo)) ||
715                        (!i_index.isEmpty() && i_index.numElements() == ni_glo*nj_glo));
716   }
717
718   // Check global zoom of a domain
719   // If there is no zoom defined for the domain, zoom will have value of global doamin
720   void CDomain::checkZoom(void)
721   {
722     if (global_zoom_ibegin.isEmpty())
723      global_zoom_ibegin.setValue(0);
724     if (global_zoom_ni.isEmpty())
725      global_zoom_ni.setValue(ni_glo);
726     if (global_zoom_jbegin.isEmpty())
727      global_zoom_jbegin.setValue(0);
728     if (global_zoom_nj.isEmpty())
729      global_zoom_nj.setValue(nj_glo);
730   }
731
732   //----------------------------------------------------------------
733
734   // Check validity of local domain on using the combination of 3 parameters: ibegin, ni and i_index
735   void CDomain::checkLocalIDomain(void)
736   {
737      // If ibegin and ni are provided then we use them to check the validity of local domain
738      if (i_index.isEmpty() && !ibegin.isEmpty() && !ni.isEmpty())
739      {
740        if ((ni.getValue() < 0 || ibegin.getValue() < 0) || ((ibegin.getValue() + ni.getValue()) > ni_glo.getValue()))
741        {
742          ERROR("CDomain::checkLocalIDomain(void)",
743                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
744                << "The local domain is wrongly defined,"
745                << " check the attributes 'ni_glo' (" << ni_glo.getValue() << "), 'ni' (" << ni.getValue() << ") and 'ibegin' (" << ibegin.getValue() << ")");
746        }
747      }
748
749      // i_index has higher priority than ibegin and ni
750      if (!i_index.isEmpty())
751      {
752        int minIIndex = (0 < i_index.numElements()) ? i_index(0) : 0;
753        if (ni.isEmpty()) 
754        {         
755         // No information about ni
756          int minIndex = ni_glo - 1;
757          int maxIndex = 0;
758          for (int idx = 0; idx < i_index.numElements(); ++idx)
759          {
760            if (i_index(idx) < minIndex) minIndex = i_index(idx);
761            if (i_index(idx) > maxIndex) maxIndex = i_index(idx);
762          }
763          ni = maxIndex - minIndex + 1; 
764          minIIndex = minIIndex;         
765        }
766
767        // It's not so correct but if ibegin is not the first value of i_index
768        // then data on local domain has user-defined distribution. In this case, ibegin, ni have no meaning.
769        if (ibegin.isEmpty()) ibegin = minIIndex;
770      }
771      else if (ibegin.isEmpty() && ni.isEmpty())
772      {
773        ibegin = 0;
774        ni = ni_glo;
775      }
776      else if ((!ibegin.isEmpty() && ni.isEmpty()) || (ibegin.isEmpty() && !ni.isEmpty()))
777      {
778        ERROR("CDomain::checkLocalIDomain(void)",
779              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
780              << "The local domain is wrongly defined," << endl
781              << "i_index is empty and either 'ni' or 'ibegin' is not defined. " 
782              << "If 'ni' and 'ibegin' are used to define a domain, both of them must not be empty.");
783      }
784       
785
786      if ((ni.getValue() < 0 || ibegin.getValue() < 0))
787      {
788        ERROR("CDomain::checkLocalIDomain(void)",
789              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
790              << "The local domain is wrongly defined,"
791              << " check the attributes 'ni_glo' (" << ni_glo.getValue() << "), 'ni' (" << ni.getValue() << ") and 'ibegin' (" << ibegin.getValue() << ")");
792      }
793   }
794
795   // Check validity of local domain on using the combination of 3 parameters: jbegin, nj and j_index
796   void CDomain::checkLocalJDomain(void)
797   {
798    // If jbegin and nj are provided then we use them to check the validity of local domain
799     if (j_index.isEmpty() && !jbegin.isEmpty() && !nj.isEmpty())
800     {
801       if ((nj.getValue() < 0 || jbegin.getValue() < 0) || (jbegin.getValue() + nj.getValue()) > nj_glo.getValue())
802       {
803         ERROR("CDomain::checkLocalJDomain(void)",
804                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
805                << "The local domain is wrongly defined,"
806                << " check the attributes 'nj_glo' (" << nj_glo.getValue() << "), 'nj' (" << nj.getValue() << ") and 'jbegin' (" << jbegin.getValue() << ")");
807       }
808     }
809
810     if (!j_index.isEmpty())
811     {
812        int minJIndex = (0 < j_index.numElements()) ? j_index(0) : 0;
813        if (nj.isEmpty()) 
814        {
815          // No information about nj
816          int minIndex = nj_glo - 1;
817          int maxIndex = 0;
818          for (int idx = 0; idx < j_index.numElements(); ++idx)
819          {
820            if (j_index(idx) < minIndex) minIndex = j_index(idx);
821            if (j_index(idx) > maxIndex) maxIndex = j_index(idx);
822          }
823          nj = maxIndex - minIndex + 1;
824          minJIndex = minIndex; 
825        } 
826        // It's the same as checkLocalIDomain. It's not so correct but if jbegin is not the first value of j_index
827        // then data on local domain has user-defined distribution. In this case, jbegin has no meaning.
828       if (jbegin.isEmpty()) jbegin = minJIndex;       
829     }
830     else if (jbegin.isEmpty() && nj.isEmpty())
831     {
832       jbegin = 0;
833       nj = nj_glo;
834     }     
835
836
837     if ((nj.getValue() < 0 || jbegin.getValue() < 0))
838     {
839       ERROR("CDomain::checkLocalJDomain(void)",
840              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
841              << "The local domain is wrongly defined,"
842              << " check the attributes 'nj_glo' (" << nj_glo.getValue() << "), 'nj' (" << nj.getValue() << ") and 'jbegin' (" << jbegin.getValue() << ")");
843     }
844   }
845
846   //----------------------------------------------------------------
847
848   void CDomain::checkMask(void)
849   {
850      if (!mask_1d.isEmpty() && !mask_2d.isEmpty())
851        ERROR("CDomain::checkMask(void)",
852              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
853              << "Both mask_1d and mask_2d are defined but only one can be used at the same time." << std::endl
854              << "Please define only one mask: 'mask_1d' or 'mask_2d'.");
855
856      if (!mask_1d.isEmpty() && mask_2d.isEmpty())
857      {
858        if (mask_1d.numElements() != i_index.numElements())
859          ERROR("CDomain::checkMask(void)",
860                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
861                << "'mask_1d' does not have the same size as the local domain." << std::endl
862                << "Local size is " << i_index.numElements() << "." << std::endl
863                << "Mask size is " << mask_1d.numElements() << ".");
864      }
865
866      if (mask_1d.isEmpty() && !mask_2d.isEmpty())
867      {
868        if (mask_2d.extent(0) != ni || mask_2d.extent(1) != nj)
869          ERROR("CDomain::checkMask(void)",
870                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
871                << "The mask does not have the same size as the local domain." << std::endl
872                << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
873                << "Mask size is " << mask_2d.extent(0) << " x " << mask_2d.extent(1) << ".");
874      }
875
876      if (!mask_2d.isEmpty())
877      {
878        mask_1d.resize(mask_2d.extent(0) * mask_2d.extent(1));
879        for (int j = 0; j < nj; ++j)
880          for (int i = 0; i < ni; ++i) mask_1d(i+j*ni) = mask_2d(i,j);
881        mask_2d.reset();
882      }
883      else if (mask_1d.isEmpty())
884      {
885        mask_1d.resize(i_index.numElements());
886        for (int i = 0; i < i_index.numElements(); ++i) mask_1d(i) = true;
887      }
888   }
889
890   //----------------------------------------------------------------
891
892   void CDomain::checkDomainData(void)
893   {
894      if (data_dim.isEmpty())
895      {
896        data_dim.setValue(1);
897      }
898      else if (!(data_dim.getValue() == 1 || data_dim.getValue() == 2))
899      {
900        ERROR("CDomain::checkDomainData(void)",
901              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
902              << "The data dimension is invalid, 'data_dim' must be 1 or 2 not << " << data_dim.getValue() << ".");
903      }
904
905      if (data_ibegin.isEmpty())
906         data_ibegin.setValue(0);
907      if (data_jbegin.isEmpty())
908         data_jbegin.setValue(0);
909
910      if (data_ni.isEmpty())
911      {
912        data_ni.setValue((data_dim == 1) ? (ni.getValue() * nj.getValue()) : ni.getValue());
913      }
914      else if (data_ni.getValue() < 0)
915      {
916        ERROR("CDomain::checkDomainData(void)",
917              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
918              << "The data size cannot be negative ('data_ni' = " << data_ni.getValue() << ").");
919      }
920
921      if (data_nj.isEmpty())
922      {
923        data_nj.setValue((data_dim.getValue() == 1) ? (ni.getValue() * nj.getValue()) : nj.getValue());
924      }
925      else if (data_nj.getValue() < 0)
926      {
927        ERROR("CDomain::checkDomainData(void)",
928              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
929              << "The data size cannot be negative ('data_nj' = " << data_nj.getValue() << ").");
930      }
931   }
932
933   //----------------------------------------------------------------
934
935   void CDomain::checkCompression(void)
936   {
937      if (!data_i_index.isEmpty())
938      {
939        if (!data_j_index.isEmpty() &&
940            data_j_index.numElements() != data_i_index.numElements())
941        {
942           ERROR("CDomain::checkCompression(void)",
943                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
944                 << "'data_i_index' and 'data_j_index' arrays must have the same size." << std::endl
945                 << "'data_i_index' size = " << data_i_index.numElements() << std::endl
946                 << "'data_j_index' size = " << data_j_index.numElements());
947        }
948
949        if (2 == data_dim)
950        {
951          if (data_j_index.isEmpty())
952          {
953             ERROR("CDomain::checkCompression(void)",
954                   << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
955                   << "'data_j_index' must be defined when 'data_i_index' is set and 'data_dim' is 2.");
956          }
957        }
958        else // (1 == data_dim)
959        {
960          if (data_j_index.isEmpty())
961          {
962            data_j_index.resize(data_ni);
963            for (int j = 0; j < data_ni; ++j) data_j_index(j) = 0;
964          }
965        }
966      }
967      else
968      {
969        if (data_dim == 2 && !data_j_index.isEmpty())
970          ERROR("CDomain::checkCompression(void)",
971                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
972                << "'data_i_index' must be defined when 'data_j_index' is set and 'data_dim' is 2.");
973
974        if (1 == data_dim)
975        {
976          data_i_index.resize(data_ni);
977          data_j_index.resize(data_ni);
978
979          for (int i = 0; i < data_ni; ++i)
980          {
981            data_i_index(i) = i;
982            data_j_index(i) = 0;
983          }
984        }
985        else // (data_dim == 2)
986        {
987          const int dsize = data_ni * data_nj;
988          data_i_index.resize(dsize);
989          data_j_index.resize(dsize);
990
991          for(int count = 0, j = 0; j < data_nj; ++j)
992          {
993            for(int i = 0; i < data_ni; ++i, ++count)
994            {
995              data_i_index(count) = i;
996              data_j_index(count) = j;
997            }
998          }
999        }
1000      }
1001   }
1002
1003   //----------------------------------------------------------------
1004   void CDomain::computeLocalMask(void)
1005   {
1006     localMask.resize(ni*nj) ;
1007     localMask=false ;
1008     size_t zoom_ibegin=global_zoom_ibegin ;
1009     size_t zoom_iend=global_zoom_ibegin+global_zoom_ni-1 ;
1010     size_t zoom_jbegin=global_zoom_jbegin ;
1011     size_t zoom_jend=global_zoom_jbegin+global_zoom_nj-1 ;
1012
1013
1014     size_t dn=data_i_index.numElements() ;
1015     int i,j ;
1016     size_t k,ind ;
1017
1018     for(k=0;k<dn;k++)
1019     {
1020       if (data_dim==2)
1021       {
1022          i=data_i_index(k)+data_ibegin ;
1023          j=data_j_index(k)+data_jbegin ;
1024       }
1025       else
1026       {
1027          i=(data_i_index(k)+data_ibegin)%ni ;
1028          j=(data_i_index(k)+data_ibegin)/ni ;
1029       }
1030
1031       if (i>=0 && i<ni && j>=0 && j<nj)
1032         if (i+ibegin>=zoom_ibegin && i+ibegin<=zoom_iend && j+jbegin>=zoom_jbegin && j+jbegin<=zoom_jend)
1033         {
1034           ind=i+ni*j ;
1035           localMask(ind)=mask_1d(ind) ;
1036         }
1037     }
1038   }
1039
1040
1041
1042
1043
1044
1045
1046   void CDomain::checkEligibilityForCompressedOutput(void)
1047   {
1048     // We don't check if the mask or the indexes are valid here, just if they have been defined at this point.
1049     isCompressible_ = !mask_1d.isEmpty() || !mask_2d.isEmpty() || !data_i_index.isEmpty();
1050   }
1051
1052   //----------------------------------------------------------------
1053
1054   void CDomain::completeLonLatClient(void)
1055   {
1056     if (!lonvalue_2d.isEmpty())
1057     {
1058       lonvalue_client.resize(ni * nj);
1059       latvalue_client.resize(ni * nj);
1060       if (hasBounds)
1061       {
1062         bounds_lon_client.resize(nvertex, ni * nj);
1063         bounds_lat_client.resize(nvertex, ni * nj);
1064       }
1065
1066       for (int j = 0; j < nj; ++j)
1067       {
1068         for (int i = 0; i < ni; ++i)
1069         {
1070           int k = j * ni + i;
1071
1072           lonvalue_client(k) = lonvalue_2d(i,j);
1073           latvalue_client(k) = latvalue_2d(i,j);
1074
1075           if (hasBounds)
1076           {
1077             for (int n = 0; n < nvertex; ++n)
1078             {
1079               bounds_lon_client(n,k) = bounds_lon_2d(n,i,j);
1080               bounds_lat_client(n,k) = bounds_lat_2d(n,i,j);
1081             }
1082           }
1083         }
1084       }
1085     }
1086     else if (!lonvalue_1d.isEmpty())
1087     {
1088       if (type_attr::rectilinear == type)
1089       {
1090         if (ni == lonvalue_1d.numElements() && nj == latvalue_1d.numElements())
1091         {
1092           lonvalue_client.resize(ni * nj);
1093           latvalue_client.resize(ni * nj);
1094           if (hasBounds)
1095           {
1096             bounds_lon_client.resize(nvertex, ni * nj);
1097             bounds_lat_client.resize(nvertex, ni * nj);
1098           }
1099
1100           for (int j = 0; j < nj; ++j)
1101           {
1102             for (int i = 0; i < ni; ++i)
1103             {
1104               int k = j * ni + i;
1105
1106               lonvalue_client(k) = lonvalue_1d(i);
1107               latvalue_client(k) = latvalue_1d(j);
1108
1109               if (hasBounds)
1110               {
1111                 for (int n = 0; n < nvertex; ++n)
1112                 {
1113                   bounds_lon_client(n,k) = bounds_lon_1d(n,i);
1114                   bounds_lat_client(n,k) = bounds_lat_1d(n,j);
1115                 }
1116               }
1117             }
1118           }
1119         }
1120         else if (i_index.numElements() == lonvalue_1d.numElements() && j_index.numElements() == latvalue_1d.numElements())
1121         {
1122           lonvalue_client.reference(lonvalue_1d);
1123           latvalue_client.reference(latvalue_1d);
1124            if (hasBounds)
1125           {
1126             bounds_lon_client.reference(bounds_lon_1d);
1127             bounds_lat_client.reference(bounds_lat_1d);
1128           }
1129         }
1130         else
1131           ERROR("CDomain::completeLonClient(void)",
1132                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1133                 << "'lonvalue_1d' and 'latvalue_1d' does not have the same size as the local domain." << std::endl
1134                 << "'lonvalue_1d' size is " << lonvalue_1d.numElements() 
1135                 << " and 'latvalue_1d' size is " << latvalue_1d.numElements() << std::endl
1136                 << " They should be correspondingly " << ni.getValue() << " and "  << nj.getValue() << " or " << std::endl
1137                 << i_index.numElements() << " and "  << j_index.numElements() << ".");
1138       }
1139       else if (type == type_attr::curvilinear || type == type_attr::unstructured)
1140       {
1141         lonvalue_client.reference(lonvalue_1d);
1142         latvalue_client.reference(latvalue_1d);
1143         if (hasBounds)
1144         {
1145           bounds_lon_client.reference(bounds_lon_1d);
1146           bounds_lat_client.reference(bounds_lat_1d);
1147         }
1148       }
1149     }
1150   }
1151
1152   void CDomain::checkBounds(void)
1153   {
1154     if (!nvertex.isEmpty() && nvertex > 0)
1155     {
1156       if (!bounds_lon_1d.isEmpty() && !bounds_lon_2d.isEmpty())
1157         ERROR("CDomain::checkBounds(void)",
1158               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1159               << "Only one longitude boundary attribute can be used but both 'bounds_lon_1d' and 'bounds_lon_2d' are defined." << std::endl
1160               << "Define only one longitude boundary attribute: 'bounds_lon_1d' or 'bounds_lon_2d'.");
1161
1162       if (!bounds_lat_1d.isEmpty() && !bounds_lat_2d.isEmpty())
1163         ERROR("CDomain::checkBounds(void)",
1164               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1165               << "Only one latitude boundary attribute can be used but both 'bounds_lat_1d' and 'bounds_lat_2d' are defined." << std::endl
1166               << "Define only one latitude boundary attribute: 'bounds_lat_1d' or 'bounds_lat_2d'.");
1167
1168       if ((!bounds_lon_1d.isEmpty() && bounds_lat_1d.isEmpty()) || (bounds_lon_1d.isEmpty() && !bounds_lat_1d.isEmpty()))
1169       {
1170         ERROR("CDomain::checkBounds(void)",
1171               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1172               << "Only 'bounds_lon_1d' or 'bounds_lat_1d' is defined." << std::endl
1173               << "Please define either both attributes or none.");
1174       }
1175
1176       if ((!bounds_lon_2d.isEmpty() && bounds_lat_2d.isEmpty()) || (bounds_lon_2d.isEmpty() && !bounds_lat_2d.isEmpty()))
1177       {
1178         ERROR("CDomain::checkBounds(void)",
1179               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1180               << "Only 'bounds_lon_2d' or 'bounds_lat_2d' is defined." << std::endl
1181               << "Please define either both attributes or none.");
1182       }
1183
1184       if (!bounds_lon_1d.isEmpty() && nvertex.getValue() != bounds_lon_1d.extent(0))
1185         ERROR("CDomain::checkBounds(void)",
1186               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1187               << "'bounds_lon_1d' dimension is not compatible with 'nvertex'." << std::endl
1188               << "'bounds_lon_1d' dimension is " << bounds_lon_1d.extent(1)
1189               << " but nvertex is " << nvertex.getValue() << ".");
1190
1191       if (!bounds_lon_2d.isEmpty() && nvertex.getValue() != bounds_lon_2d.extent(0))
1192         ERROR("CDomain::checkBounds(void)",
1193               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1194               << "'bounds_lon_2d' dimension is not compatible with 'nvertex'." << std::endl
1195               << "'bounds_lon_2d' dimension is " << bounds_lon_2d.extent(2)
1196               << " but nvertex is " << nvertex.getValue() << ".");
1197
1198       if (!bounds_lon_1d.isEmpty() && lonvalue_1d.isEmpty())
1199         ERROR("CDomain::checkBounds(void)",
1200               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1201               << "Since 'bounds_lon_1d' is defined, 'lonvalue_1d' must be defined too." << std::endl);
1202
1203       if (!bounds_lon_2d.isEmpty() && lonvalue_2d.isEmpty())
1204         ERROR("CDomain::checkBounds(void)",
1205               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1206               << "Since 'bounds_lon_2d' is defined, 'lonvalue_2d' must be defined too." << std::endl);
1207
1208       if (!bounds_lat_1d.isEmpty() && nvertex.getValue() != bounds_lat_1d.extent(0))
1209         ERROR("CDomain::checkBounds(void)",
1210               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1211               << "'bounds_lat_1d' dimension is not compatible with 'nvertex'." << std::endl
1212               << "'bounds_lat_1d' dimension is " << bounds_lat_1d.extent(1)
1213               << " but nvertex is " << nvertex.getValue() << ".");
1214
1215       if (!bounds_lat_2d.isEmpty() && nvertex.getValue() != bounds_lat_2d.extent(0))
1216         ERROR("CDomain::checkBounds(void)",
1217               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1218               << "'bounds_lat_2d' dimension is not compatible with 'nvertex'." << std::endl
1219               << "'bounds_lat_2d' dimension is " << bounds_lat_2d.extent(2)
1220               << " but nvertex is " << nvertex.getValue() << ".");
1221
1222       if (!bounds_lat_1d.isEmpty() && latvalue_1d.isEmpty())
1223         ERROR("CDomain::checkBounds(void)",
1224               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1225               << "Since 'bounds_lat_1d' is defined, 'latvalue_1d' must be defined too." << std::endl);
1226
1227       if (!bounds_lat_2d.isEmpty() && latvalue_2d.isEmpty())
1228         ERROR("CDomain::checkBounds(void)",
1229               << "Since 'bounds_lat_2d' is defined, 'latvalue_2d' must be defined too." << std::endl);
1230
1231       hasBounds = true;
1232     }
1233     else
1234     {
1235       hasBounds = false;
1236       nvertex = 0;
1237     }
1238   }
1239
1240   void CDomain::checkArea(void)
1241   {
1242     hasArea = !area.isEmpty();
1243     if (hasArea)
1244     {
1245       if (area.extent(0) != ni || area.extent(1) != nj)
1246       {
1247         ERROR("CDomain::checkArea(void)",
1248               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1249               << "The area does not have the same size as the local domain." << std::endl
1250               << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
1251               << "Area size is " << area.extent(0) << " x " << area.extent(1) << ".");
1252       }
1253     }
1254   }
1255
1256   void CDomain::checkLonLat()
1257   {
1258     hasLonLat = (!latvalue_1d.isEmpty() && !lonvalue_1d.isEmpty()) ||
1259                 (!latvalue_2d.isEmpty() && !lonvalue_2d.isEmpty());
1260     if (hasLonLat)
1261     {
1262       if (!lonvalue_1d.isEmpty() && !lonvalue_2d.isEmpty())
1263         ERROR("CDomain::checkLonLat()",
1264               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1265               << "Only one longitude attribute can be used but both 'lonvalue_1d' and 'lonvalue_2d' are defined." << std::endl
1266               << "Define only one longitude attribute: 'lonvalue_1d' or 'lonvalue_2d'.");
1267
1268       if (!lonvalue_1d.isEmpty() && lonvalue_2d.isEmpty())
1269       {
1270         if ((type_attr::rectilinear != type) && (lonvalue_1d.numElements() != i_index.numElements()))
1271           ERROR("CDomain::checkLonLat()",
1272                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1273                 << "'lonvalue_1d' does not have the same size as the local domain." << std::endl
1274                 << "Local size is " << i_index.numElements() << "." << std::endl
1275                 << "'lonvalue_1d' size is " << lonvalue_1d.numElements() << ".");
1276       }
1277
1278       if (lonvalue_1d.isEmpty() && !lonvalue_2d.isEmpty())
1279       {
1280         if (lonvalue_2d.extent(0) != ni || lonvalue_2d.extent(1) != nj)
1281           ERROR("CDomain::checkLonLat()",
1282                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1283                 << "'lonvalue_2d' does not have the same size as the local domain." << std::endl
1284                 << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
1285                 << "'lonvalue_2d' size is " << lonvalue_2d.extent(0) << " x " << lonvalue_2d.extent(1) << ".");
1286       }
1287
1288       if (!latvalue_1d.isEmpty() && !latvalue_2d.isEmpty())
1289         ERROR("CDomain::checkLonLat()",
1290               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1291               << "Only one latitude attribute can be used but both 'latvalue_1d' and 'latvalue_2d' are defined." << std::endl
1292               << "Define only one latitude attribute: 'latvalue_1d' or 'latvalue_2d'.");
1293
1294       if (!latvalue_1d.isEmpty() && latvalue_2d.isEmpty())
1295       {
1296         if ((type_attr::rectilinear != type) && (latvalue_1d.numElements() != i_index.numElements()))
1297           ERROR("CDomain::checkLonLat()",
1298                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1299                 << "'latvalue_1d' does not have the same size as the local domain." << std::endl
1300                 << "Local size is " << i_index.numElements() << "." << std::endl
1301                 << "'latvalue_1d' size is " << latvalue_1d.numElements() << ".");
1302       }
1303
1304       if (latvalue_1d.isEmpty() && !latvalue_2d.isEmpty())
1305       {
1306         if (latvalue_2d.extent(0) != ni || latvalue_2d.extent(1) != nj)
1307           ERROR("CDomain::checkLonLat()",
1308                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1309                 << "'latvalue_2d' does not have the same size as the local domain." << std::endl
1310                 << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
1311                 << "'latvalue_2d' size is " << latvalue_2d.extent(0) << " x " << latvalue_2d.extent(1) << ".");
1312       }
1313     }
1314   }
1315
1316   void CDomain::checkAttributesOnClientAfterTransformation()
1317   {
1318     
1319     CContext* context=CContext::getCurrent() ;
1320
1321     if (this->isClientAfterTransformationChecked) return;
1322     if (context->hasClient)
1323     {
1324       this->checkMask(); 
1325       if (hasLonLat || hasArea || isCompressible_) 
1326        {
1327          this->computeConnectedServer();
1328        }
1329       if (hasLonLat) 
1330        {
1331          this->completeLonLatClient();
1332        }
1333     }
1334
1335     this->isClientAfterTransformationChecked = true;
1336   }
1337
1338   //----------------------------------------------------------------
1339   // Divide function checkAttributes into 2 seperate ones
1340   // This function only checks all attributes of current domain
1341   void CDomain::checkAttributesOnClient()
1342   {
1343     if (this->isClientChecked) return;
1344     CContext* context=CContext::getCurrent();
1345
1346      this->checkDomain();
1347      this->checkBounds();
1348      this->checkArea();
1349      this->checkLonLat();
1350
1351      if (context->hasClient)
1352      { // CÃŽté client uniquement
1353         this->checkMask();
1354         this->checkDomainData();
1355         this->checkCompression();
1356         this->computeLocalMask() ;
1357      }
1358      else
1359      { // CÃŽté serveur uniquement
1360      }
1361
1362      this->isClientChecked = true;
1363   }
1364
1365   // Send all checked attributes to server
1366   void CDomain::sendCheckedAttributes()
1367   {
1368     if (!this->isClientChecked) checkAttributesOnClient();
1369     if (!this->isClientAfterTransformationChecked) checkAttributesOnClientAfterTransformation();
1370     CContext* context=CContext::getCurrent() ;
1371
1372     if (this->isChecked) return;
1373     if (context->hasClient)
1374     {
1375       sendServerAttribut();
1376       if (hasLonLat || hasArea || isCompressible_) sendLonLatArea();
1377     }
1378     this->isChecked = true;
1379   }
1380
1381   void CDomain::checkAttributes(void)
1382   {
1383      if (this->isChecked) return;
1384      CContext* context=CContext::getCurrent() ;
1385
1386      this->checkDomain();
1387      this->checkLonLat();
1388      this->checkBounds();
1389      this->checkArea();
1390
1391      if (context->hasClient)
1392      { // CÃŽté client uniquement
1393         this->checkMask();
1394         this->checkDomainData();
1395         this->checkCompression();
1396         this->computeLocalMask() ;
1397
1398      }
1399      else
1400      { // CÃŽté serveur uniquement
1401      }
1402
1403      if (context->hasClient)
1404      {
1405        this->computeConnectedServer();
1406        this->completeLonLatClient();
1407        this->sendServerAttribut();
1408        this->sendLonLatArea();
1409      }
1410
1411      this->isChecked = true;
1412   }
1413
1414  void CDomain::sendServerAttribut(void)
1415  {
1416    CContext* context = CContext::getCurrent();
1417    CContextClient* client = context->client;
1418    int nbServer = client->serverSize;
1419
1420    CServerDistributionDescription serverDescription(nGlobDomain_, nbServer);
1421    if (isUnstructed_) serverDescription.computeServerDistribution(false, 0);
1422    else serverDescription.computeServerDistribution(false, 1);
1423
1424    std::vector<std::vector<int> > serverIndexBegin = serverDescription.getServerIndexBegin();
1425    std::vector<std::vector<int> > serverDimensionSizes = serverDescription.getServerDimensionSizes();
1426
1427    CEventClient event(getType(),EVENT_ID_SERVER_ATTRIBUT);
1428    if (client->isServerLeader())
1429    {
1430      std::list<CMessage> msgs;
1431
1432      const std::list<int>& ranks = client->getRanksServerLeader();
1433      for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1434      {
1435        // Use const int to ensure CMessage holds a copy of the value instead of just a reference
1436        const int ibegin_srv = serverIndexBegin[*itRank][0];
1437        const int jbegin_srv = serverIndexBegin[*itRank][1];
1438        const int ni_srv = serverDimensionSizes[*itRank][0];
1439        const int nj_srv = serverDimensionSizes[*itRank][1];
1440        const int iend_srv = ibegin_srv + ni_srv - 1;
1441        const int jend_srv = jbegin_srv + nj_srv - 1;
1442
1443        msgs.push_back(CMessage());
1444        CMessage& msg = msgs.back();
1445        msg << this->getId() ;
1446        msg << ni_srv << ibegin_srv << iend_srv << nj_srv << jbegin_srv << jend_srv;
1447        msg << global_zoom_ni.getValue() << global_zoom_ibegin.getValue() << global_zoom_nj.getValue() << global_zoom_jbegin.getValue();
1448        msg << isCompressible_;
1449
1450        event.push(*itRank,1,msg);
1451      }
1452      client->sendEvent(event);
1453    }
1454    else client->sendEvent(event);
1455  }
1456
1457  void CDomain::computeNGlobDomain()
1458  {
1459    nGlobDomain_.resize(2);
1460    nGlobDomain_[0] = ni_glo.getValue();
1461    nGlobDomain_[1] = nj_glo.getValue();
1462  }
1463
1464  void CDomain::computeConnectedServer(void)
1465  {
1466
1467    CContext* context=CContext::getCurrent() ;
1468    CContextClient* client=context->client ;
1469    int nbServer=client->serverSize;
1470    int rank = client->clientRank;
1471    bool doComputeGlobalIndexServer = true;
1472
1473    int i,j,i_ind,j_ind, nbIndex;
1474    int global_zoom_iend=global_zoom_ibegin+global_zoom_ni-1 ;
1475    int global_zoom_jend=global_zoom_jbegin+global_zoom_nj-1 ;
1476
1477    // Precompute number of index
1478    int globalIndexCountZoom = 0;
1479    nbIndex = i_index.numElements();
1480    for (i = 0; i < nbIndex; ++i)
1481    {
1482      i_ind=i_index(i);
1483      j_ind=j_index(i);
1484
1485      if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1486      {
1487        ++globalIndexCountZoom;
1488      }
1489    }
1490
1491    int globalIndexWrittenCount = 0;
1492    if (isCompressible_)
1493    {
1494      for (i = 0; i < data_i_index.numElements(); ++i)
1495      {
1496        i_ind = CDistributionClient::getDomainIndex(data_i_index(i), data_j_index(i),
1497                                                    data_ibegin, data_jbegin, data_dim, ni,
1498                                                    j_ind);
1499        if (i_ind >= 0 && i_ind < ni && j_ind >= 0 && j_ind < nj && mask_1d(i_ind + j_ind * ni))
1500        {
1501          i_ind += ibegin;
1502          j_ind += jbegin;
1503          if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1504            ++globalIndexWrittenCount;
1505        }
1506      }
1507    }
1508
1509    // Fill in index
1510    CArray<size_t,1> globalIndexDomainZoom(globalIndexCountZoom);
1511    CArray<size_t,1> localIndexDomainZoom(globalIndexCountZoom);
1512    CArray<size_t,1> globalIndexDomain(nbIndex);
1513    size_t globalIndex;
1514    int globalIndexCount = 0;
1515    globalIndexCountZoom = 0;
1516
1517    for (i = 0; i < nbIndex; ++i)
1518    {
1519      i_ind=i_index(i);
1520      j_ind=j_index(i);
1521      globalIndex = i_ind + j_ind * ni_glo;
1522      globalIndexDomain(globalIndexCount) = globalIndex;
1523      ++globalIndexCount;
1524      if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1525      {
1526        globalIndexDomainZoom(globalIndexCountZoom) = globalIndex;
1527        localIndexDomainZoom(globalIndexCountZoom) = i;
1528        ++globalIndexCountZoom;
1529      }
1530    }
1531
1532    CArray<int,1> globalIndexWrittenDomain(globalIndexWrittenCount);
1533    if (isCompressible_)
1534    {
1535      globalIndexWrittenCount = 0;
1536      for (i = 0; i < data_i_index.numElements(); ++i)
1537      {
1538        i_ind = CDistributionClient::getDomainIndex(data_i_index(i), data_j_index(i),
1539                                                    data_ibegin, data_jbegin, data_dim, ni,
1540                                                    j_ind);
1541        if (i_ind >= 0 && i_ind < ni && j_ind >= 0 && j_ind < nj && mask_1d(i_ind + j_ind * ni))
1542        {
1543          i_ind += ibegin;
1544          j_ind += jbegin;
1545          if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1546          {
1547            globalIndexWrittenDomain(globalIndexWrittenCount) = i_ind + j_ind * ni_glo;
1548            ++globalIndexWrittenCount;
1549          }
1550        }
1551      }
1552    }
1553   
1554
1555    size_t globalSizeIndex = 1, indexBegin, indexEnd;
1556    int range, clientSize = client->clientSize;
1557    for (int i = 0; i < nGlobDomain_.size(); ++i) globalSizeIndex *= nGlobDomain_[i];
1558    indexBegin = 0;
1559    if (globalSizeIndex <= clientSize)
1560    {
1561      indexBegin = rank%globalSizeIndex;
1562      indexEnd = indexBegin;
1563    }
1564    else
1565    {
1566      for (int i = 0; i < clientSize; ++i)
1567      {
1568        range = globalSizeIndex / clientSize;
1569        if (i < (globalSizeIndex%clientSize)) ++range;
1570        if (i == client->clientRank) break;
1571        indexBegin += range;
1572      }
1573      indexEnd = indexBegin + range - 1;
1574    }
1575   
1576    CServerDistributionDescription serverDescription(nGlobDomain_, nbServer);
1577    if (isUnstructed_) serverDescription.computeServerGlobalIndexInRange(std::make_pair<size_t,size_t>(indexBegin, indexEnd), 0);
1578    else serverDescription.computeServerGlobalIndexInRange(std::make_pair<size_t,size_t>(indexBegin, indexEnd), 1);
1579
1580    CClientServerMapping* clientServerMap = new CClientServerMappingDistributed(serverDescription.getGlobalIndexRange(),
1581                                                                                client->intraComm);
1582    clientServerMap->computeServerIndexMapping(globalIndexDomain); 
1583   
1584    const CClientServerMapping::GlobalIndexMap& globalIndexDomainOnServer = clientServerMap->getGlobalIndexOnServer();
1585   
1586    CClientServerMapping::GlobalIndexMap::const_iterator it  = globalIndexDomainOnServer.begin(),
1587                                                         ite = globalIndexDomainOnServer.end();
1588    typedef XIOSBinarySearchWithIndex<size_t> BinarySearch;
1589    std::vector<int>::iterator itVec;
1590
1591    indSrv_.clear();
1592    indWrittenSrv_.clear();
1593    for (; it != ite; ++it)
1594    {
1595      int rank = it->first;
1596      int indexSize = it->second.size();
1597      std::vector<int> permutIndex(indexSize);
1598      XIOSAlgorithms::fillInIndex(indexSize, permutIndex);
1599      XIOSAlgorithms::sortWithIndex<size_t, CVectorStorage>(it->second, permutIndex);
1600      BinarySearch binSearch(it->second);
1601      int nb = globalIndexDomainZoom.numElements();
1602      for (int i = 0; i < nb; ++i)
1603      {
1604        if (binSearch.search(permutIndex.begin(), permutIndex.end(), globalIndexDomainZoom(i), itVec))
1605        {
1606          indSrv_[rank].push_back(localIndexDomainZoom(i));
1607        }
1608      }
1609      for (int i = 0; i < globalIndexWrittenDomain.numElements(); ++i)
1610      {
1611        if (binSearch.search(permutIndex.begin(), permutIndex.end(), globalIndexWrittenDomain(i), itVec))
1612        {
1613          indWrittenSrv_[rank].push_back(globalIndexWrittenDomain(i));
1614        }
1615      }
1616    }
1617
1618    connectedServerRank_.clear();
1619    for (it = globalIndexDomainOnServer.begin(); it != ite; ++it) {
1620      connectedServerRank_.push_back(it->first);
1621    }
1622
1623    nbConnectedClients_ = clientServerMap->computeConnectedClients(client->serverSize, client->clientSize, client->intraComm, connectedServerRank_);
1624
1625    delete clientServerMap;
1626  }
1627
1628  const std::map<int, vector<size_t> >& CDomain::getIndexServer() const
1629  {
1630    return indSrv_;
1631  }
1632
1633  /*!
1634    Send index from client to server(s)
1635  */
1636  void CDomain::sendIndex()
1637  {
1638    int ns, n, i, j, ind, nv, idx;
1639    CContext* context = CContext::getCurrent();
1640    CContextClient* client=context->client;
1641
1642    CEventClient eventIndex(getType(), EVENT_ID_INDEX);
1643
1644    list<CMessage> list_msgsIndex;
1645    list<CArray<int,1> > list_indi, list_indj, list_writtenInd;
1646
1647    std::map<int, std::vector<size_t> >::const_iterator it, iteMap;
1648    iteMap = indSrv_.end();
1649    for (int k = 0; k < connectedServerRank_.size(); ++k)
1650    {
1651      int nbData = 0;
1652      int rank = connectedServerRank_[k];
1653      it = indSrv_.find(rank);
1654      if (iteMap != it)
1655        nbData = it->second.size();
1656
1657      list_indi.push_back(CArray<int,1>(nbData));
1658      list_indj.push_back(CArray<int,1>(nbData));
1659
1660      CArray<int,1>& indi = list_indi.back();
1661      CArray<int,1>& indj = list_indj.back();
1662      const std::vector<size_t>& temp = it->second;
1663      for (n = 0; n < nbData; ++n)
1664      {
1665        idx = static_cast<int>(it->second[n]);
1666        indi(n) = i_index(idx);
1667        indj(n) = j_index(idx);
1668      }
1669
1670      list_msgsIndex.push_back(CMessage());
1671
1672      list_msgsIndex.back() << this->getId() << (int)type; // enum ne fonctionne pour les message => ToFix
1673      list_msgsIndex.back() << isCurvilinear;
1674      list_msgsIndex.back() << list_indi.back() << list_indj.back();
1675
1676      if (isCompressible_)
1677      {
1678        std::vector<int>& writtenIndSrc = indWrittenSrv_[rank];
1679        list_writtenInd.push_back(CArray<int,1>(writtenIndSrc.size()));
1680        CArray<int,1>& writtenInd = list_writtenInd.back();
1681
1682        for (n = 0; n < writtenInd.numElements(); ++n)
1683          writtenInd(n) = writtenIndSrc[n];
1684
1685        list_msgsIndex.back() << writtenInd;
1686      }
1687
1688      eventIndex.push(rank, nbConnectedClients_[rank], list_msgsIndex.back());
1689    }
1690
1691    client->sendEvent(eventIndex);
1692  }
1693
1694  /*!
1695    Send area from client to server(s)
1696  */
1697  void CDomain::sendArea()
1698  {
1699    if (!hasArea) return;
1700
1701    int ns, n, i, j, ind, nv, idx;
1702    CContext* context = CContext::getCurrent();
1703    CContextClient* client=context->client;
1704
1705    // send area for each connected server
1706    CEventClient eventArea(getType(), EVENT_ID_AREA);
1707
1708    list<CMessage> list_msgsArea;
1709    list<CArray<double,1> > list_area;
1710
1711    std::map<int, std::vector<size_t> >::const_iterator it, iteMap;
1712    iteMap = indSrv_.end();
1713    for (int k = 0; k < connectedServerRank_.size(); ++k)
1714    {
1715      int nbData = 0;
1716      int rank = connectedServerRank_[k];
1717      it = indSrv_.find(rank);
1718      if (iteMap != it)
1719        nbData = it->second.size();
1720      list_area.push_back(CArray<double,1>(nbData));
1721
1722      const std::vector<size_t>& temp = it->second;
1723      for (n = 0; n < nbData; ++n)
1724      {
1725        idx = static_cast<int>(it->second[n]);
1726        i = i_index(idx);
1727        j = j_index(idx);
1728        if (hasArea)
1729          list_area.back()(n) = area(i - ibegin, j - jbegin);
1730      }
1731
1732      list_msgsArea.push_back(CMessage());
1733      list_msgsArea.back() << this->getId() << list_area.back();
1734      eventArea.push(rank, nbConnectedClients_[rank], list_msgsArea.back());
1735    }
1736    client->sendEvent(eventArea);
1737  }
1738
1739  /*!
1740    Send longitude and latitude from client to servers
1741    Each client send long and lat information to corresponding connected server(s).
1742    Because longitude and latitude are optional, this function only called if latitude and longitude exist
1743  */
1744  void CDomain::sendLonLat()
1745  {
1746    if (!hasLonLat) return;
1747
1748    int ns, n, i, j, ind, nv, idx;
1749    CContext* context = CContext::getCurrent();
1750    CContextClient* client=context->client;
1751
1752    // send lon lat for each connected server
1753    CEventClient eventLon(getType(), EVENT_ID_LON);
1754    CEventClient eventLat(getType(), EVENT_ID_LAT);
1755
1756    list<CMessage> list_msgsLon, list_msgsLat;
1757    list<CArray<double,1> > list_lon, list_lat;
1758    list<CArray<double,2> > list_boundslon, list_boundslat;
1759
1760    std::map<int, std::vector<size_t> >::const_iterator it, iteMap;
1761    iteMap = indSrv_.end();
1762    for (int k = 0; k < connectedServerRank_.size(); ++k)
1763    {
1764      int nbData = 0;
1765      int rank = connectedServerRank_[k];
1766      it = indSrv_.find(rank);
1767      if (iteMap != it)
1768        nbData = it->second.size();
1769
1770      list_lon.push_back(CArray<double,1>(nbData));
1771      list_lat.push_back(CArray<double,1>(nbData));
1772
1773      if (hasBounds)
1774      {
1775        list_boundslon.push_back(CArray<double,2>(nvertex, nbData));
1776        list_boundslat.push_back(CArray<double,2>(nvertex, nbData));
1777      }
1778
1779      CArray<double,1>& lon = list_lon.back();
1780      CArray<double,1>& lat = list_lat.back();
1781      const std::vector<size_t>& temp = it->second;
1782      for (n = 0; n < nbData; ++n)
1783      {
1784        idx = static_cast<int>(it->second[n]);
1785        lon(n) = lonvalue_client(idx);
1786        lat(n) = latvalue_client(idx);
1787
1788        if (hasBounds)
1789        {
1790          CArray<double,2>& boundslon = list_boundslon.back();
1791          CArray<double,2>& boundslat = list_boundslat.back();
1792
1793          for (nv = 0; nv < nvertex; ++nv)
1794          {
1795            boundslon(nv, n) = bounds_lon_client(nv, idx);
1796            boundslat(nv, n) = bounds_lat_client(nv, idx);
1797          }
1798        }
1799      }
1800
1801      list_msgsLon.push_back(CMessage());
1802      list_msgsLat.push_back(CMessage());
1803
1804      list_msgsLon.back() << this->getId() << list_lon.back();
1805      list_msgsLat.back() << this->getId() << list_lat.back();
1806
1807      if (hasBounds)
1808      {
1809        list_msgsLon.back() << list_boundslon.back();
1810        list_msgsLat.back() << list_boundslat.back();
1811      }
1812
1813      eventLon.push(rank, nbConnectedClients_[rank], list_msgsLon.back());
1814      eventLat.push(rank, nbConnectedClients_[rank], list_msgsLat.back());
1815    }
1816
1817    client->sendEvent(eventLon);
1818    client->sendEvent(eventLat);
1819  }
1820
1821  /*!
1822    Send some optional information to server(s)
1823    In the future, this function can be extended with more optional information to send
1824  */
1825  void CDomain::sendLonLatArea(void)
1826  {
1827    sendIndex();
1828    sendLonLat();
1829    sendArea();
1830  }
1831
1832  bool CDomain::dispatchEvent(CEventServer& event)
1833  {
1834    if (SuperClass::dispatchEvent(event)) return true;
1835    else
1836    {
1837      switch(event.type)
1838      {
1839        case EVENT_ID_SERVER_ATTRIBUT:
1840          recvServerAttribut(event);
1841          return true;
1842          break;
1843        case EVENT_ID_INDEX:
1844          recvIndex(event);
1845          return true;
1846          break;
1847        case EVENT_ID_LON:
1848          recvLon(event);
1849          return true;
1850          break;
1851        case EVENT_ID_LAT:
1852          recvLat(event);
1853          return true;
1854          break;
1855        case EVENT_ID_AREA:
1856          recvArea(event);
1857          return true;
1858          break;
1859        default:
1860          ERROR("bool CDomain::dispatchEvent(CEventServer& event)",
1861                << "Unknown Event");
1862          return false;
1863       }
1864    }
1865  }
1866
1867  /*!
1868    Receive attributes event from clients(s)
1869    \param[in] event event contain info about rank and associated attributes
1870  */
1871  void CDomain::recvServerAttribut(CEventServer& event)
1872  {
1873    CBufferIn* buffer=event.subEvents.begin()->buffer;
1874    string domainId ;
1875    *buffer>>domainId ;
1876    get(domainId)->recvServerAttribut(*buffer) ;
1877  }
1878
1879  /*!
1880    Receive attributes from client(s): zoom info and begin and n of each server
1881    \param[in] rank rank of client source
1882    \param[in] buffer message containing attributes info
1883  */
1884  void CDomain::recvServerAttribut(CBufferIn& buffer)
1885  {
1886    int global_zoom_ni_tmp, global_zoom_ibegin_tmp, global_zoom_nj_tmp, global_zoom_jbegin_tmp;
1887    buffer >> ni_srv >> ibegin_srv >> iend_srv >> nj_srv >> jbegin_srv >> jend_srv
1888           >> global_zoom_ni_tmp >> global_zoom_ibegin_tmp >> global_zoom_nj_tmp >> global_zoom_jbegin_tmp
1889           >> isCompressible_;
1890
1891    global_zoom_ni.setValue(global_zoom_ni_tmp);
1892    global_zoom_ibegin.setValue(global_zoom_ibegin_tmp);
1893    global_zoom_nj.setValue(global_zoom_nj_tmp);
1894    global_zoom_jbegin.setValue(global_zoom_jbegin_tmp);
1895
1896    int zoom_iend = global_zoom_ibegin + global_zoom_ni - 1;
1897    int zoom_jend = global_zoom_jbegin + global_zoom_nj - 1;
1898
1899    zoom_ibegin_srv = global_zoom_ibegin > ibegin_srv ? global_zoom_ibegin : ibegin_srv ;
1900    zoom_iend_srv = zoom_iend < iend_srv ? zoom_iend : iend_srv ;
1901    zoom_ni_srv=zoom_iend_srv-zoom_ibegin_srv+1 ;
1902
1903    zoom_jbegin_srv = global_zoom_jbegin > jbegin_srv ? global_zoom_jbegin : jbegin_srv ;
1904    zoom_jend_srv = zoom_jend < jend_srv ? zoom_jend : jend_srv ;
1905    zoom_nj_srv=zoom_jend_srv-zoom_jbegin_srv+1 ;
1906
1907    if (zoom_ni_srv<=0 || zoom_nj_srv<=0)
1908    {
1909      zoom_ibegin_srv=0 ; zoom_iend_srv=0 ; zoom_ni_srv=0 ;
1910      zoom_jbegin_srv=0 ; zoom_jend_srv=0 ; zoom_nj_srv=0 ;
1911    }
1912    lonvalue_srv.resize(zoom_ni_srv*zoom_nj_srv) ;
1913    lonvalue_srv = 0. ;
1914    latvalue_srv.resize(zoom_ni_srv*zoom_nj_srv) ;
1915    latvalue_srv = 0. ;
1916    if (hasBounds)
1917    {
1918      bounds_lon_srv.resize(nvertex,zoom_ni_srv*zoom_nj_srv) ;
1919      bounds_lon_srv = 0. ;
1920      bounds_lat_srv.resize(nvertex,zoom_ni_srv*zoom_nj_srv) ;
1921      bounds_lat_srv = 0. ;
1922    }
1923
1924    if (hasArea)
1925    {
1926      area_srv.resize(zoom_ni_srv * zoom_nj_srv);
1927      area_srv = 0.;
1928    }
1929
1930  }
1931
1932  /*!
1933    Receive index event from clients(s)
1934    \param[in] event event contain info about rank and associated index
1935  */
1936  void CDomain::recvIndex(CEventServer& event)
1937  {
1938    CDomain* domain;
1939
1940    list<CEventServer::SSubEvent>::iterator it;
1941    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
1942    {
1943      CBufferIn* buffer = it->buffer;
1944      string domainId;
1945      *buffer >> domainId;
1946      domain = get(domainId);
1947      domain->recvIndex(it->rank, *buffer);
1948    }
1949
1950    if (domain->isCompressible_)
1951    {
1952      std::sort(domain->indexesToWrite.begin(), domain->indexesToWrite.end());
1953
1954      CContextServer* server = CContext::getCurrent()->server;
1955      domain->numberWrittenIndexes_ = domain->indexesToWrite.size();
1956      MPI_Allreduce(&domain->numberWrittenIndexes_, &domain->totalNumberWrittenIndexes_, 1, MPI_INT, MPI_SUM, server->intraComm);
1957      MPI_Scan(&domain->numberWrittenIndexes_, &domain->offsetWrittenIndexes_, 1, MPI_INT, MPI_SUM, server->intraComm);
1958      domain->offsetWrittenIndexes_ -= domain->numberWrittenIndexes_;
1959    }
1960  }
1961
1962  /*!
1963    Receive index information from client(s)
1964    \param[in] rank rank of client source
1965    \param[in] buffer message containing index info
1966  */
1967  void CDomain::recvIndex(int rank, CBufferIn& buffer)
1968  {
1969    int type_int;
1970    buffer >> type_int >> isCurvilinear >> indiSrv[rank] >> indjSrv[rank];
1971    type.setValue((type_attr::t_enum)type_int); // probleme des type enum avec les buffers : ToFix
1972
1973    if (isCompressible_)
1974    {
1975      CArray<int, 1> writtenIndexes;
1976      buffer >> writtenIndexes;
1977      indexesToWrite.reserve(indexesToWrite.size() + writtenIndexes.numElements());
1978      for (int i = 0; i < writtenIndexes.numElements(); ++i)
1979        indexesToWrite.push_back(writtenIndexes(i));
1980    }
1981  }
1982
1983  /*!
1984    Receive longitude event from clients(s)
1985    \param[in] event event contain info about rank and associated longitude
1986  */
1987  void CDomain::recvLon(CEventServer& event)
1988  {
1989    list<CEventServer::SSubEvent>::iterator it;
1990    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
1991    {
1992      CBufferIn* buffer = it->buffer;
1993      string domainId;
1994      *buffer >> domainId;
1995      get(domainId)->recvLon(it->rank, *buffer);
1996    }
1997  }
1998
1999  /*!
2000    Receive longitude information from client(s)
2001    \param[in] rank rank of client source
2002    \param[in] buffer message containing longitude info
2003  */
2004  void CDomain::recvLon(int rank, CBufferIn& buffer)
2005  {
2006    CArray<int,1> &indi = indiSrv[rank], &indj = indjSrv[rank];
2007    CArray<double,1> lon;
2008    CArray<double,2> boundslon;
2009
2010    buffer >> lon;
2011
2012    if (hasBounds) buffer >> boundslon;
2013
2014    int i, j, ind_srv;
2015    for (int ind = 0; ind < indi.numElements(); ind++)
2016    {
2017      i = indi(ind); j = indj(ind);
2018      ind_srv = (i - zoom_ibegin_srv) + (j - zoom_jbegin_srv) * zoom_ni_srv;
2019      lonvalue_srv(ind_srv) = lon(ind);
2020      if (hasBounds)
2021      {
2022        for (int nv = 0; nv < nvertex; ++nv)
2023          bounds_lon_srv(nv, ind_srv) = boundslon(nv, ind);
2024      }
2025    }
2026  }
2027
2028  /*!
2029    Receive latitude event from clients(s)
2030    \param[in] event event contain info about rank and associated latitude
2031  */
2032  void CDomain::recvLat(CEventServer& event)
2033  {
2034    list<CEventServer::SSubEvent>::iterator it;
2035    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2036    {
2037      CBufferIn* buffer = it->buffer;
2038      string domainId;
2039      *buffer >> domainId;
2040      get(domainId)->recvLat(it->rank, *buffer);
2041    }
2042  }
2043
2044  /*!
2045    Receive latitude information from client(s)
2046    \param[in] rank rank of client source
2047    \param[in] buffer message containing latitude info
2048  */
2049  void CDomain::recvLat(int rank, CBufferIn& buffer)
2050  {
2051    CArray<int,1> &indi = indiSrv[rank], &indj = indjSrv[rank];
2052    CArray<double,1> lat;
2053    CArray<double,2> boundslat;
2054
2055    buffer >> lat;
2056    if (hasBounds) buffer >> boundslat;
2057
2058    int i, j, ind_srv;
2059    for (int ind = 0; ind < indi.numElements(); ind++)
2060    {
2061      i = indi(ind); j = indj(ind);
2062      ind_srv = (i - zoom_ibegin_srv) + (j - zoom_jbegin_srv) * zoom_ni_srv;
2063      latvalue_srv(ind_srv) = lat(ind);
2064      if (hasBounds)
2065      {
2066        for (int nv = 0; nv < nvertex; nv++)
2067          bounds_lat_srv(nv, ind_srv) = boundslat(nv, ind);
2068      }
2069    }
2070  }
2071
2072  /*!
2073    Receive area event from clients(s)
2074    \param[in] event event contain info about rank and associated area
2075  */
2076  void CDomain::recvArea(CEventServer& event)
2077  {
2078    list<CEventServer::SSubEvent>::iterator it;
2079    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2080    {
2081      CBufferIn* buffer = it->buffer;
2082      string domainId;
2083      *buffer >> domainId;
2084      get(domainId)->recvArea(it->rank, *buffer);
2085    }
2086  }
2087
2088  /*!
2089    Receive area information from client(s)
2090    \param[in] rank rank of client source
2091    \param[in] buffer message containing area info
2092  */
2093  void CDomain::recvArea(int rank, CBufferIn& buffer)
2094  {
2095    CArray<int,1> &indi = indiSrv[rank], &indj = indjSrv[rank];
2096    CArray<double,1> clientArea;
2097
2098    buffer >> clientArea;
2099
2100    int i, j, ind_srv;
2101    for (int ind = 0; ind < indi.numElements(); ind++)
2102    {
2103      i = indi(ind); j = indj(ind);
2104      ind_srv = (i - zoom_ibegin_srv) + (j - zoom_jbegin_srv) * zoom_ni_srv;
2105      area_srv(ind_srv) = clientArea(ind);
2106    }
2107  }
2108
2109  CTransformation<CDomain>* CDomain::addTransformation(ETranformationType transType, const StdString& id)
2110  {
2111    transformationMap_.push_back(std::make_pair(transType, CTransformation<CDomain>::createTransformation(transType,id)));
2112    return transformationMap_.back().second;
2113  }
2114
2115  /*!
2116    Check whether a domain has transformation
2117    \return true if domain has transformation
2118  */
2119  bool CDomain::hasTransformation()
2120  {
2121    return (!transformationMap_.empty());
2122  }
2123
2124  /*!
2125    Set transformation for current domain. It's the method to move transformation in hierarchy
2126    \param [in] domTrans transformation on domain
2127  */
2128  void CDomain::setTransformations(const TransMapTypes& domTrans)
2129  {
2130    transformationMap_ = domTrans;
2131  }
2132
2133  /*!
2134    Get all transformation current domain has
2135    \return all transformation
2136  */
2137  CDomain::TransMapTypes CDomain::getAllTransformations(void)
2138  {
2139    return transformationMap_;
2140  }
2141
2142  /*!
2143    Check the validity of all transformations applied on domain
2144  This functions is called AFTER all inherited attributes are solved
2145  */
2146  void CDomain::checkTransformations()
2147  {
2148    TransMapTypes::const_iterator itb = transformationMap_.begin(), it,
2149                                  ite = transformationMap_.end();
2150//    for (it = itb; it != ite; ++it)
2151//    {
2152//      (it->second)->checkValid(this);
2153//    }
2154  }
2155
2156  void CDomain::duplicateTransformation(CDomain* src)
2157  {
2158    if (src->hasTransformation())
2159    {
2160      this->setTransformations(src->getAllTransformations());
2161    }
2162  }
2163
2164  /*!
2165   * Go through the hierarchy to find the domain from which the transformations must be inherited
2166   */
2167  void CDomain::solveInheritanceTransformation()
2168  {
2169    if (hasTransformation() || !hasDirectDomainReference())
2170      return;
2171
2172    CDomain* domain = this;
2173    std::vector<CDomain*> refDomains;
2174    while (!domain->hasTransformation() && domain->hasDirectDomainReference())
2175    {
2176      refDomains.push_back(domain);
2177      domain = domain->getDirectDomainReference();
2178    }
2179
2180    if (domain->hasTransformation())
2181      for (size_t i = 0; i < refDomains.size(); ++i)
2182        refDomains[i]->setTransformations(domain->getAllTransformations());
2183  }
2184
2185  /*!
2186    Parse children nodes of a domain in xml file.
2187    Whenver there is a new transformation, its type and name should be added into this function
2188    \param node child node to process
2189  */
2190  void CDomain::parse(xml::CXMLNode & node)
2191  {
2192    SuperClass::parse(node);
2193
2194    if (node.goToChildElement())
2195    {
2196      StdString nodeElementName;
2197      do
2198      {
2199        StdString nodeId("");
2200        if (node.getAttributes().end() != node.getAttributes().find("id"))
2201        { nodeId = node.getAttributes()["id"]; }
2202
2203        nodeElementName = node.getElementName();
2204        if(transformationMapList_ptr == 0) initializeTransformationMap();
2205        std::map<StdString, ETranformationType>::const_iterator ite = (*transformationMapList_ptr).end(), it;
2206        it = (*transformationMapList_ptr).find(nodeElementName);
2207        if (ite != it)
2208        {
2209          transformationMap_.push_back(std::make_pair(it->second, CTransformation<CDomain>::createTransformation(it->second,
2210                                                                                                                nodeId,
2211                                                                                                                &node)));
2212        }
2213        else
2214        {
2215          ERROR("void CDomain::parse(xml::CXMLNode & node)",
2216                << "The transformation " << nodeElementName << " has not been supported yet.");
2217        }
2218      } while (node.goToNextElement()) ;
2219      node.goToParentElement();
2220    }
2221  }
2222   //----------------------------------------------------------------
2223
2224   DEFINE_REF_FUNC(Domain,domain)
2225
2226   ///---------------------------------------------------------------
2227
2228} // namespace xios
Note: See TracBrowser for help on using the repository browser.