source: XIOS/dev/dev_olga/src/node/domain.cpp @ 1132

Last change on this file since 1132 was 1132, checked in by mhnguyen, 7 years ago

Correcting a minor bug on writting unstructured grid

+) Correct the mapping between received data and written data.
+) Clean some redundant codes

Test
+) On Curie
+) Writing on unstructured grid works correctly

  • 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: 112.6 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   /// ////////////////////// Dfinitions ////////////////////// ///
29
30   CDomain::CDomain(void)
31      : CObjectTemplate<CDomain>(), CDomainAttributes()
32      , isChecked(false), relFiles(), isClientChecked(false), nbConnectedClients_(), indSrv_(), connectedServerRank_(), connectedServerZoomRank_()
33      , hasBounds(false), hasArea(false), isDistributed_(false), isCompressible_(false), isUnstructed_(false)
34      , isClientAfterTransformationChecked(false), hasLonLat(false)
35      , isRedistributed_(false), hasPole(false), doZoomByIndex_(false)
36      , lonvalue(), latvalue(), bounds_lonvalue(), bounds_latvalue()
37      , globalLocalIndexMap_(), computedWrittenIndex_(false)
38   {
39   }
40
41   CDomain::CDomain(const StdString & id)
42      : CObjectTemplate<CDomain>(id), CDomainAttributes()
43      , isChecked(false), relFiles(), isClientChecked(false), nbConnectedClients_(), indSrv_(), connectedServerRank_(), connectedServerZoomRank_()
44      , hasBounds(false), hasArea(false), isDistributed_(false), isCompressible_(false), isUnstructed_(false)
45      , isClientAfterTransformationChecked(false), hasLonLat(false)
46      , isRedistributed_(false), hasPole(false), doZoomByIndex_(false)
47      , lonvalue(), latvalue(), bounds_lonvalue(), bounds_latvalue()
48      , globalLocalIndexMap_(), computedWrittenIndex_(false)
49   {
50         }
51
52   CDomain::~CDomain(void)
53   {
54   }
55
56   ///---------------------------------------------------------------
57
58   void CDomain::assignMesh(const StdString meshName, const int nvertex)
59   {
60     mesh = CMesh::getMesh(meshName, nvertex);
61   }
62
63   CDomain* CDomain::createDomain()
64   {
65     CDomain* domain = CDomainGroup::get("domain_definition")->createChild();
66     return domain;
67   }
68
69   std::map<StdString, ETranformationType> CDomain::transformationMapList_ = std::map<StdString, ETranformationType>();
70   bool CDomain::_dummyTransformationMapList = CDomain::initializeTransformationMap(CDomain::transformationMapList_);
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   const std::set<StdString> & CDomain::getRelFiles(void) const
82   {
83      return (this->relFiles);
84   }
85
86
87   const std::vector<int>& CDomain::getIndexesToWrite(void) const
88   {
89     return indexesToWrite;
90   }
91
92   /*!
93     Returns the number of indexes written by each server.
94     \return the number of indexes written by each server
95   */
96   int CDomain::getNumberWrittenIndexes() const
97   {
98     return numberWrittenIndexes_;
99   }
100
101   /*!
102     Returns the total number of indexes written by the servers.
103     \return the total number of indexes written by the servers
104   */
105   int CDomain::getTotalNumberWrittenIndexes() const
106   {
107     return totalNumberWrittenIndexes_;
108   }
109
110   /*!
111     Returns the offset of indexes written by each server.
112     \return the offset of indexes written by each server
113   */
114   int CDomain::getOffsetWrittenIndexes() const
115   {
116     return offsetWrittenIndexes_;
117   }
118
119   /*!
120     Returns the start of indexes written by each server.
121     \return the start of indexes written by each server
122   */
123   const std::vector<int>& CDomain::getStartWriteIndex() const
124   {
125     return start_write_index_;
126   }
127
128   /*!
129     Returns the count of indexes written by each server.
130     \return the count of indexes written by each server
131   */
132   const std::vector<int>& CDomain::getCountWriteIndex() const
133   {
134     return count_write_index_;
135   }
136
137   /*!
138     Returns the local data written by each server.     
139   */
140   const std::vector<int>& CDomain::getLocalWriteSize() const
141   {
142     return local_write_size_;
143   }
144
145   /*!
146     Returns the global data written by all server.     
147   */
148   const std::vector<int>& CDomain::getGlobalWriteSize() const
149   {
150     return global_write_size_;
151   }
152
153   //----------------------------------------------------------------
154
155   /*!
156    * Compute the minimum buffer size required to send the attributes to the server(s).
157    *
158    * \return A map associating the server rank with its minimum buffer size.
159    */
160   std::map<int, StdSize> CDomain::getAttributesBufferSize()
161   {
162     CContext* context = CContext::getCurrent();
163     // For now the assumption is that secondary server pools consist of the same number of procs.
164     // CHANGE the line below if the assumption changes.
165     CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[0] : context->client;
166
167     std::map<int, StdSize> attributesSizes = getMinimumBufferSizeForAttributes();
168
169     if (client->isServerLeader())
170     {
171       // size estimation for sendDistributionAttribut
172       size_t size = 11 * sizeof(size_t);
173
174       const std::list<int>& ranks = client->getRanksServerLeader();
175       for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
176       {
177         if (size > attributesSizes[*itRank])
178           attributesSizes[*itRank] = size;
179       }
180     }
181
182     boost::unordered_map<int, vector<size_t> >::const_iterator itIndexEnd = indSrv_.end();
183     std::map<int, std::vector<int> >::const_iterator itWrittenIndexEnd = indWrittenSrv_.end();
184     for (size_t k = 0; k < connectedServerRank_.size(); ++k)
185     {
186       int rank = connectedServerRank_[k];
187       boost::unordered_map<int, std::vector<size_t> >::const_iterator it = indSrv_.find(rank);
188       size_t idxCount = (it != itIndexEnd) ? it->second.size() : 0;
189
190       // size estimation for sendIndex (and sendArea which is always smaller or equal)
191       size_t sizeIndexEvent = 2 * sizeof(size_t) + 2 * CArray<int,1>::size(idxCount);
192       if (isCompressible_)
193       {
194         std::map<int, std::vector<int> >::const_iterator itWritten = indWrittenSrv_.find(rank);
195         size_t writtenIdxCount = (itWritten != itWrittenIndexEnd) ? itWritten->second.size() : 0;
196         sizeIndexEvent += CArray<int,1>::size(writtenIdxCount);
197       }
198
199       // size estimation for sendLonLat
200       size_t sizeLonLatEvent = CArray<double,1>::size(idxCount);
201       if (hasBounds)
202         sizeLonLatEvent += CArray<double,2>::size(nvertex * idxCount);
203
204       size_t size = CEventClient::headerSize + getId().size() + sizeof(size_t) + std::max(sizeIndexEvent, sizeLonLatEvent);
205       if (size > attributesSizes[rank])
206         attributesSizes[rank] = size;
207     }
208
209     return attributesSizes;
210   }
211
212   //----------------------------------------------------------------
213
214   bool CDomain::isEmpty(void) const
215   {
216      return ((this->zoom_i_index.isEmpty()) || (0 == this->zoom_i_index.numElements()));
217      // return ((this->zoom_ni_srv == 0) ||
218      //         (this->zoom_nj_srv == 0));
219   }
220
221   //----------------------------------------------------------------
222
223   bool CDomain::IsWritten(const StdString & filename) const
224   {
225      return (this->relFiles.find(filename) != this->relFiles.end());
226   }
227
228   bool CDomain::isWrittenCompressed(const StdString& filename) const
229   {
230      return (this->relFilesCompressed.find(filename) != this->relFilesCompressed.end());
231   }
232
233   //----------------------------------------------------------------
234
235   bool CDomain::isDistributed(void) const
236   {
237      return isDistributed_;
238   }
239
240   //----------------------------------------------------------------
241
242   /*!
243    * Test whether the data defined on the domain can be outputted in a compressed way.
244    *
245    * \return true if and only if a mask was defined for this domain
246    */
247   bool CDomain::isCompressible(void) const
248   {
249      return isCompressible_;
250   }
251
252   void CDomain::addRelFile(const StdString & filename)
253   {
254      this->relFiles.insert(filename);
255   }
256
257   void CDomain::addRelFileCompressed(const StdString& filename)
258   {
259      this->relFilesCompressed.insert(filename);
260   }
261
262   StdString CDomain::GetName(void)   { return (StdString("domain")); }
263   StdString CDomain::GetDefName(void){ return (CDomain::GetName()); }
264   ENodeType CDomain::GetType(void)   { return (eDomain); }
265
266   //----------------------------------------------------------------
267
268   /*!
269     Redistribute RECTILINEAR domain with a number of local domains.
270   All attributes ni,nj,ibegin,jbegin (if defined) will be rewritten
271   The optional attributes lonvalue, latvalue will be added. Because this function only serves (for now)
272   for interpolation from unstructured domain to rectilinear one, range of latvalue is 0-360 and lonvalue is -90 - +90
273    \param [in] nbLocalDomain number of local domain on the domain destination
274   */
275   void CDomain::redistribute(int nbLocalDomain)
276   {
277     if (this->isRedistributed_) return;
278
279     this->isRedistributed_ = true;
280     CContext* context = CContext::getCurrent();
281     // For now the assumption is that secondary server pools consist of the same number of procs.
282     // CHANGE the line below if the assumption changes.
283     CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[0] : context->client;
284     int rankClient = client->clientRank;
285     int rankOnDomain = rankClient%nbLocalDomain;
286
287     if (ni_glo.isEmpty() || ni_glo <= 0 )
288     {
289        ERROR("CDomain::redistribute(int nbLocalDomain)",
290           << "[ Id = " << this->getId() << " ] "
291           << "The global domain is badly defined,"
292           << " check the \'ni_glo\'  value !")
293     }
294
295     if (nj_glo.isEmpty() || nj_glo <= 0 )
296     {
297        ERROR("CDomain::redistribute(int nbLocalDomain)",
298           << "[ Id = " << this->getId() << " ] "
299           << "The global domain is badly defined,"
300           << " check the \'nj_glo\'  value !")
301     }
302
303     if ((type_attr::rectilinear == type)  || (type_attr::curvilinear == type))
304     {
305        int globalDomainSize = ni_glo * nj_glo;
306        if (globalDomainSize <= nbLocalDomain)
307        {
308          for (int idx = 0; idx < nbLocalDomain; ++idx)
309          {
310            if (rankOnDomain < globalDomainSize)
311            {
312              int iIdx = rankOnDomain % ni_glo;
313              int jIdx = rankOnDomain / ni_glo;
314              ibegin.setValue(iIdx); jbegin.setValue(jIdx);
315              ni.setValue(1); nj.setValue(1);
316            }
317            else
318            {
319              ibegin.setValue(0); jbegin.setValue(0);
320              ni.setValue(0); nj.setValue(0);
321            }
322          }
323        }
324        else
325        {
326          float njGlo = nj_glo.getValue();
327          float niGlo = ni_glo.getValue();
328          int nbProcOnX, nbProcOnY, range;
329
330          // Compute (approximately) number of segment on x and y axis
331          float yOverXRatio = njGlo/niGlo;
332
333          nbProcOnX = std::ceil(std::sqrt(nbLocalDomain/yOverXRatio));
334          nbProcOnY = std::ceil(((float)nbLocalDomain)/nbProcOnX);
335
336          // Simple distribution: Sweep from top to bottom, left to right
337          // Calculate local begin on x
338          std::vector<int> ibeginVec(nbProcOnX,0), jbeginVec(nbProcOnY,0);
339          std::vector<int> niVec(nbProcOnX), njVec(nbProcOnY);
340          for (int i = 1; i < nbProcOnX; ++i)
341          {
342            range = ni_glo / nbProcOnX;
343            if (i < (ni_glo%nbProcOnX)) ++range;
344            niVec[i-1] = range;
345            ibeginVec[i] = ibeginVec[i-1] + niVec[i-1];
346          }
347          niVec[nbProcOnX-1] = ni_glo - ibeginVec[nbProcOnX-1];
348
349          // Calculate local begin on y
350          for (int j = 1; j < nbProcOnY; ++j)
351          {
352            range = nj_glo / nbProcOnY;
353            if (j < (nj_glo%nbProcOnY)) ++range;
354            njVec[j-1] = range;
355            jbeginVec[j] = jbeginVec[j-1] + njVec[j-1];
356          }
357          njVec[nbProcOnY-1] = nj_glo - jbeginVec[nbProcOnY-1];
358
359          // Now assign value to ni, ibegin, nj, jbegin
360          int iIdx = rankOnDomain % nbProcOnX;
361          int jIdx = rankOnDomain / nbProcOnX;
362
363          if (rankOnDomain != (nbLocalDomain-1))
364          {
365            ibegin.setValue(ibeginVec[iIdx]);
366            jbegin.setValue(jbeginVec[jIdx]);
367            nj.setValue(njVec[jIdx]);
368            ni.setValue(niVec[iIdx]);
369          }
370          else // just merge all the remaining rectangle into the last one
371          {
372            ibegin.setValue(ibeginVec[iIdx]);
373            jbegin.setValue(jbeginVec[jIdx]);
374            nj.setValue(njVec[jIdx]);
375            ni.setValue(ni_glo - ibeginVec[iIdx]);
376          }
377        }
378
379        // Now fill other attributes
380        if (type_attr::rectilinear == type) fillInRectilinearLonLat();
381     }
382     else  // unstructured domain
383     {
384       if (this->i_index.isEmpty())
385       {
386          int globalDomainSize = ni_glo * nj_glo;
387          if (globalDomainSize <= nbLocalDomain)
388          {
389            for (int idx = 0; idx < nbLocalDomain; ++idx)
390            {
391              if (rankOnDomain < globalDomainSize)
392              {
393                int iIdx = rankOnDomain % ni_glo;
394                int jIdx = rankOnDomain / ni_glo;
395                ibegin.setValue(iIdx); jbegin.setValue(jIdx);
396                ni.setValue(1); nj.setValue(1);
397              }
398              else
399              {
400                ibegin.setValue(0); jbegin.setValue(0);
401                ni.setValue(0); nj.setValue(0);
402              }
403            }
404          }
405          else
406          {
407            float njGlo = nj_glo.getValue();
408            float niGlo = ni_glo.getValue();
409            std::vector<int> ibeginVec(nbLocalDomain,0);
410            std::vector<int> niVec(nbLocalDomain);
411            for (int i = 1; i < nbLocalDomain; ++i)
412            {
413              int range = ni_glo / nbLocalDomain;
414              if (i < (ni_glo%nbLocalDomain)) ++range;
415              niVec[i-1] = range;
416              ibeginVec[i] = ibeginVec[i-1] + niVec[i-1];
417            }
418            niVec[nbLocalDomain-1] = ni_glo - ibeginVec[nbLocalDomain-1];
419
420            int iIdx = rankOnDomain % nbLocalDomain;
421            ibegin.setValue(ibeginVec[iIdx]);
422            jbegin.setValue(0);
423            ni.setValue(niVec[iIdx]);
424            nj.setValue(1);
425          }
426        }
427        else
428        {
429          ibegin.setValue(this->i_index(0));
430          jbegin.setValue(0);
431          ni.setValue(this->i_index.numElements());
432          nj.setValue(1);
433        }
434     }
435
436     checkDomain();
437   }
438
439   /*!
440     Fill in the values for lonvalue_1d and latvalue_1d of rectilinear domain
441     Range of longitude value from 0 - 360
442     Range of latitude value from -90 - +90
443   */
444   void CDomain::fillInRectilinearLonLat()
445   {
446     if (!lonvalue_rectilinear_read_from_file.isEmpty())
447     {
448       lonvalue_1d.resize(ni);
449       for (int idx = 0; idx < ni; ++idx)
450         lonvalue_1d(idx) = lonvalue_rectilinear_read_from_file(idx+ibegin);
451       lon_start.setValue(lonvalue_rectilinear_read_from_file(0));
452       lon_end.setValue(lonvalue_rectilinear_read_from_file(ni_glo-1));
453     }
454     else
455     {
456       if (!lonvalue_2d.isEmpty()) lonvalue_2d.free();
457       lonvalue_1d.resize(ni);
458       double lonRange = lon_end - lon_start;
459       double lonStep = (1 == ni_glo.getValue()) ? lonRange : lonRange/double(ni_glo.getValue()-1);
460
461        // Assign lon value
462       for (int i = 0; i < ni; ++i)
463       {
464         if (0 == (ibegin + i))
465         {
466           lonvalue_1d(i) = lon_start;
467         }
468         else if (ni_glo == (ibegin + i + 1))
469         {
470           lonvalue_1d(i) = lon_end;
471         }
472         else
473         {
474           lonvalue_1d(i) = (ibegin + i) * lonStep  + lon_start;
475         }
476       }
477     }
478
479
480     if (!latvalue_rectilinear_read_from_file.isEmpty())
481     {
482       latvalue_1d.resize(nj);
483       for (int idx = 0; idx < nj; ++idx)
484         latvalue_1d(idx) = latvalue_rectilinear_read_from_file(idx+jbegin);
485       lat_start.setValue(latvalue_rectilinear_read_from_file(0));
486       lat_end.setValue(latvalue_rectilinear_read_from_file(nj_glo-1));
487     }
488     else
489     {
490       if (!latvalue_2d.isEmpty()) latvalue_1d.free();
491       latvalue_1d.resize(nj);
492
493       double latRange = lat_end - lat_start;
494       double latStep = (1 == nj_glo.getValue()) ? latRange : latRange/double(nj_glo.getValue()-1);
495
496       for (int j = 0; j < nj; ++j)
497       {
498         if (0 == (jbegin + j))
499         {
500            latvalue_1d(j) = lat_start;
501         }
502         else if (nj_glo == (jbegin + j + 1))
503         {
504            latvalue_1d(j) = lat_end;
505         }
506         else
507         {
508           latvalue_1d(j) =  (jbegin + j) * latStep + lat_start;
509         }
510       }
511     }
512   }
513
514
515
516   void CDomain::AllgatherRectilinearLonLat(CArray<double,1>& lon, CArray<double,1>& lat, CArray<double,1>& lon_g, CArray<double,1>& lat_g)
517   {
518          CContext* context = CContext::getCurrent();
519    // For now the assumption is that secondary server pools consist of the same number of procs.
520    // CHANGE the line below if the assumption changes.
521    CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[0] : context->client;
522          lon_g.resize(ni_glo) ;
523          lat_g.resize(nj_glo) ;
524
525
526          int* ibegin_g = new int[client->clientSize] ;
527          int* jbegin_g = new int[client->clientSize] ;
528          int* ni_g = new int[client->clientSize] ;
529          int* nj_g = new int[client->clientSize] ;
530          int v ;
531          v=ibegin ;
532          MPI_Allgather(&v,1,MPI_INT,ibegin_g,1,MPI_INT,client->intraComm) ;
533          v=jbegin ;
534          MPI_Allgather(&v,1,MPI_INT,jbegin_g,1,MPI_INT,client->intraComm) ;
535          v=ni ;
536          MPI_Allgather(&v,1,MPI_INT,ni_g,1,MPI_INT,client->intraComm) ;
537          v=nj ;
538          MPI_Allgather(&v,1,MPI_INT,nj_g,1,MPI_INT,client->intraComm) ;
539
540          MPI_Allgatherv(lon.dataFirst(),ni,MPI_DOUBLE,lon_g.dataFirst(),ni_g, ibegin_g,MPI_DOUBLE,client->intraComm) ;
541          MPI_Allgatherv(lat.dataFirst(),nj,MPI_DOUBLE,lat_g.dataFirst(),nj_g, jbegin_g,MPI_DOUBLE,client->intraComm) ;
542
543      delete[] ibegin_g ;
544      delete[] jbegin_g ;
545      delete[] ni_g ;
546      delete[] nj_g ;
547   }
548
549   void CDomain::fillInRectilinearBoundLonLat(CArray<double,1>& lon, CArray<double,1>& lat,
550                                              CArray<double,2>& boundsLon, CArray<double,2>& boundsLat)
551   {
552     int i,j,k;
553
554     const int nvertexValue = 4;
555     boundsLon.resize(nvertexValue,ni*nj);
556
557     if (ni_glo>1)
558     {
559       double lonStepStart = lon(1)-lon(0);
560       bounds_lon_start=lon(0) - lonStepStart/2;
561       double lonStepEnd = lon(ni_glo-1)-lon(ni_glo-2);
562       bounds_lon_end=lon(ni_glo-1) + lonStepEnd/2;
563       double errorBoundsLon = std::abs(360-std::abs(bounds_lon_end-bounds_lon_start));
564
565       // if errorBoundsLon is reasonably small (0.1 x cell size) consider it as closed in longitude
566       if (errorBoundsLon < std::abs(lonStepStart)*1e-1 || errorBoundsLon < std::abs(lonStepEnd)*1e-1 )
567       {
568         bounds_lon_start= (lon(0) + lon(ni_glo-1)-360)/2 ;
569         bounds_lon_end= (lon(0) +360 + lon(ni_glo-1))/2 ;
570       }
571     }
572     else
573     {
574       if (bounds_lon_start.isEmpty()) bounds_lon_start=-180. ;
575       if (bounds_lon_end.isEmpty()) bounds_lon_end=180.-1e-8 ;
576     }
577
578     for(j=0;j<nj;++j)
579       for(i=0;i<ni;++i)
580       {
581         k=j*ni+i;
582         boundsLon(0,k) = boundsLon(1,k) = (0 == (ibegin + i)) ? bounds_lon_start
583                                                               : (lon(ibegin + i)+lon(ibegin + i-1))/2;
584         boundsLon(2,k) = boundsLon(3,k) = ((ibegin + i + 1) == ni_glo) ? bounds_lon_end
585                                                                        : (lon(ibegin + i + 1)+lon(ibegin + i))/2;
586       }
587
588
589    boundsLat.resize(nvertexValue,nj*ni);
590    bool isNorthPole=false ;
591    bool isSouthPole=false ;
592    if (std::abs(90 - std::abs(lat(0))) < NumTraits<double>::epsilon()) isNorthPole = true;
593    if (std::abs(90 - std::abs(lat(nj_glo-1))) < NumTraits<double>::epsilon()) isSouthPole = true;
594
595    // lat boundaries beyond pole the assimilate it to pole
596    // lat boundarie is relativelly close to pole (0.1 x cell size) assimilate it to pole
597    if (nj_glo>1)
598    {
599      double latStepStart = lat(1)-lat(0);
600      if (isNorthPole) bounds_lat_start=lat(0);
601      else
602      {
603        bounds_lat_start=lat(0)-latStepStart/2;
604        if (bounds_lat_start >= 90 ) bounds_lat_start=90 ;
605        else if (bounds_lat_start <= -90 ) bounds_lat_start=-90 ;
606        else if (bounds_lat_start <= 90 && bounds_lat_start >= lat(0))
607        {
608          if ( std::abs(90-bounds_lat_start) <= 0.1*std::abs(latStepStart)) bounds_lat_start=90 ;
609        }
610        else if (bounds_lat_start >= -90 && bounds_lat_start <= lat(0))
611        {
612          if ( std::abs(-90 - bounds_lat_start) <= 0.1*std::abs(latStepStart)) bounds_lat_start=-90 ;
613        }
614      }
615
616      double latStepEnd = lat(nj_glo-1)-lat(nj_glo-2);
617      if (isSouthPole) bounds_lat_end=lat(nj_glo-1);
618      else
619      {
620        bounds_lat_end=lat(nj_glo-1)+latStepEnd/2;
621
622        if (bounds_lat_end >= 90 ) bounds_lat_end=90 ;
623        else if (bounds_lat_end <= -90 ) bounds_lat_end=-90 ;
624        else if (bounds_lat_end <= 90 && bounds_lat_end >= lat(nj_glo-1))
625        {
626          if ( std::abs(90-bounds_lat_end) <= 0.1*std::abs(latStepEnd)) bounds_lat_end=90 ;
627        }
628        else if (bounds_lat_end >= -90 && bounds_lat_end <= lat(nj_glo-1))
629        {
630          if ( std::abs(-90 - bounds_lat_end) <= 0.1*std::abs(latStepEnd)) bounds_lat_end=-90 ;
631        }
632      }
633    }
634    else
635    {
636      if (bounds_lat_start.isEmpty()) bounds_lon_start=-90. ;
637      if (bounds_lat_end.isEmpty()) bounds_lat_end=90 ;
638    }
639
640    for(j=0;j<nj;++j)
641      for(i=0;i<ni;++i)
642      {
643        k=j*ni+i;
644        boundsLat(1,k) = boundsLat(2,k) = (0 == (jbegin + j)) ? bounds_lat_start
645                                                              : (lat(jbegin + j)+lat(jbegin + j-1))/2;
646        boundsLat(0,k) = boundsLat(3,k) = ((jbegin + j +1) == nj_glo) ? bounds_lat_end
647                                                                      : (lat(jbegin + j + 1)+lat(jbegin + j))/2;
648      }
649   }
650
651   void CDomain::checkDomain(void)
652   {
653     if (type.isEmpty())
654     {
655       ERROR("CDomain::checkDomain(void)",
656             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
657             << "The domain type is mandatory, "
658             << "please define the 'type' attribute.")
659     }
660
661     if (type == type_attr::gaussian) 
662     {
663           hasPole=true ;
664             type.setValue(type_attr::unstructured) ;
665           }
666           else if (type == type_attr::rectilinear) hasPole=true ;
667         
668     if (type == type_attr::unstructured)
669     {
670        if (ni_glo.isEmpty())
671        {
672          ERROR("CDomain::checkDomain(void)",
673                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
674                << "The global domain is badly defined, "
675                << "the mandatory 'ni_glo' attribute is missing.")
676        }
677        else if (ni_glo <= 0)
678        {
679          ERROR("CDomain::checkDomain(void)",
680                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
681                << "The global domain is badly defined, "
682                << "'ni_glo' attribute should be strictly positive so 'ni_glo = " << ni_glo.getValue() << "' is invalid.")
683        }
684        isUnstructed_ = true;
685        nj_glo = 1;
686        nj = 1;
687        jbegin = 0;
688        if (!i_index.isEmpty()) ni = i_index.numElements();
689        j_index.resize(ni);
690        for(int i=0;i<ni;++i) j_index(i)=0;
691
692        if (!area.isEmpty())
693          area.transposeSelf(1, 0);
694     }
695
696     if (ni_glo.isEmpty())
697     {
698       ERROR("CDomain::checkDomain(void)",
699             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
700             << "The global domain is badly defined, "
701             << "the mandatory 'ni_glo' attribute is missing.")
702     }
703     else if (ni_glo <= 0)
704     {
705       ERROR("CDomain::checkDomain(void)",
706             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
707             << "The global domain is badly defined, "
708             << "'ni_glo' attribute should be strictly positive so 'ni_glo = " << ni_glo.getValue() << "' is invalid.")
709     }
710
711     if (nj_glo.isEmpty())
712     {
713       ERROR("CDomain::checkDomain(void)",
714             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
715             << "The global domain is badly defined, "
716             << "the mandatory 'nj_glo' attribute is missing.")
717     }
718     else if (nj_glo <= 0)
719     {
720       ERROR("CDomain::checkDomain(void)",
721             << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
722             << "The global domain is badly defined, "
723             << "'nj_glo' attribute should be strictly positive so 'nj_glo = " << nj_glo.getValue() << "' is invalid.")
724     }
725
726     checkLocalIDomain();
727     checkLocalJDomain();
728
729     if (i_index.isEmpty())
730     {
731       i_index.resize(ni*nj);
732       for (int j = 0; j < nj; ++j)
733         for (int i = 0; i < ni; ++i) i_index(i+j*ni) = i+ibegin;
734     }
735
736     if (j_index.isEmpty())
737     {
738       j_index.resize(ni*nj);
739       for (int j = 0; j < nj; ++j)
740         for (int i = 0; i < ni; ++i) j_index(i+j*ni) = j+jbegin;
741     }     
742     checkZoom();
743
744     isDistributed_ = !((!ni.isEmpty() && (ni == ni_glo) && !nj.isEmpty() && (nj == nj_glo)) ||
745                        (!i_index.isEmpty() && i_index.numElements() == ni_glo*nj_glo));
746   }
747
748   // Check global zoom of a domain
749   // If there is no zoom defined for the domain, zoom will have value of global doamin
750   void CDomain::checkZoom(void)
751   {
752     if (global_zoom_ibegin.isEmpty())
753      global_zoom_ibegin.setValue(0);
754     if (global_zoom_ni.isEmpty())
755      global_zoom_ni.setValue(ni_glo);
756     if (global_zoom_jbegin.isEmpty())
757      global_zoom_jbegin.setValue(0);
758     if (global_zoom_nj.isEmpty())
759      global_zoom_nj.setValue(nj_glo);
760    if (zoom_i_index.isEmpty()) zoom_i_index.setValue(i_index.getValue());
761    if (zoom_j_index.isEmpty()) zoom_j_index.setValue(j_index.getValue());
762    if (zoom_ibegin.isEmpty()) zoom_ibegin.setValue(ibegin);
763    if (zoom_ni.isEmpty()) zoom_ni.setValue(ni);
764    if (zoom_jbegin.isEmpty()) zoom_jbegin.setValue(jbegin);
765    if (zoom_nj.isEmpty()) zoom_nj.setValue(nj);
766   }
767
768   //----------------------------------------------------------------
769
770   // Check validity of local domain on using the combination of 3 parameters: ibegin, ni and i_index
771   void CDomain::checkLocalIDomain(void)
772   {
773      // If ibegin and ni are provided then we use them to check the validity of local domain
774      if (i_index.isEmpty() && !ibegin.isEmpty() && !ni.isEmpty())
775      {
776        if ((ni.getValue() < 0 || ibegin.getValue() < 0) || ((ibegin.getValue() + ni.getValue()) > ni_glo.getValue()))
777        {
778          ERROR("CDomain::checkLocalIDomain(void)",
779                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
780                << "The local domain is wrongly defined,"
781                << " check the attributes 'ni_glo' (" << ni_glo.getValue() << "), 'ni' (" << ni.getValue() << ") and 'ibegin' (" << ibegin.getValue() << ")");
782        }
783      }
784
785      // i_index has higher priority than ibegin and ni
786      if (!i_index.isEmpty())
787      {
788        int minIIndex = (0 < i_index.numElements()) ? i_index(0) : 0;
789        if (ni.isEmpty()) 
790        {         
791         // No information about ni
792          int minIndex = ni_glo - 1;
793          int maxIndex = 0;
794          for (int idx = 0; idx < i_index.numElements(); ++idx)
795          {
796            if (i_index(idx) < minIndex) minIndex = i_index(idx);
797            if (i_index(idx) > maxIndex) maxIndex = i_index(idx);
798          }
799          ni = maxIndex - minIndex + 1; 
800          minIIndex = minIIndex;         
801        }
802
803        // It's not so correct but if ibegin is not the first value of i_index
804        // then data on local domain has user-defined distribution. In this case, ibegin, ni have no meaning.
805        if (ibegin.isEmpty()) ibegin = minIIndex;
806      }
807      else if (ibegin.isEmpty() && ni.isEmpty())
808      {
809        ibegin = 0;
810        ni = ni_glo;
811      }
812      else if ((!ibegin.isEmpty() && ni.isEmpty()) || (ibegin.isEmpty() && !ni.isEmpty()))
813      {
814        ERROR("CDomain::checkLocalIDomain(void)",
815              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
816              << "The local domain is wrongly defined," << endl
817              << "i_index is empty and either 'ni' or 'ibegin' is not defined. " 
818              << "If 'ni' and 'ibegin' are used to define a domain, both of them must not be empty.");
819      }
820       
821
822      if ((ni.getValue() < 0 || ibegin.getValue() < 0))
823      {
824        ERROR("CDomain::checkLocalIDomain(void)",
825              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
826              << "The local domain is wrongly defined,"
827              << " check the attributes 'ni_glo' (" << ni_glo.getValue() << "), 'ni' (" << ni.getValue() << ") and 'ibegin' (" << ibegin.getValue() << ")");
828      }
829   }
830
831   // Check validity of local domain on using the combination of 3 parameters: jbegin, nj and j_index
832   void CDomain::checkLocalJDomain(void)
833   {
834    // If jbegin and nj are provided then we use them to check the validity of local domain
835     if (j_index.isEmpty() && !jbegin.isEmpty() && !nj.isEmpty())
836     {
837       if ((nj.getValue() < 0 || jbegin.getValue() < 0) || (jbegin.getValue() + nj.getValue()) > nj_glo.getValue())
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     if (!j_index.isEmpty())
847     {
848        int minJIndex = (0 < j_index.numElements()) ? j_index(0) : 0;
849        if (nj.isEmpty()) 
850        {
851          // No information about nj
852          int minIndex = nj_glo - 1;
853          int maxIndex = 0;
854          for (int idx = 0; idx < j_index.numElements(); ++idx)
855          {
856            if (j_index(idx) < minIndex) minIndex = j_index(idx);
857            if (j_index(idx) > maxIndex) maxIndex = j_index(idx);
858          }
859          nj = maxIndex - minIndex + 1;
860          minJIndex = minIndex; 
861        } 
862        // It's the same as checkLocalIDomain. It's not so correct but if jbegin is not the first value of j_index
863        // then data on local domain has user-defined distribution. In this case, jbegin has no meaning.
864       if (jbegin.isEmpty()) jbegin = minJIndex;       
865     }
866     else if (jbegin.isEmpty() && nj.isEmpty())
867     {
868       jbegin = 0;
869       nj = nj_glo;
870     }     
871
872
873     if ((nj.getValue() < 0 || jbegin.getValue() < 0))
874     {
875       ERROR("CDomain::checkLocalJDomain(void)",
876              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
877              << "The local domain is wrongly defined,"
878              << " check the attributes 'nj_glo' (" << nj_glo.getValue() << "), 'nj' (" << nj.getValue() << ") and 'jbegin' (" << jbegin.getValue() << ")");
879     }
880   }
881
882   //----------------------------------------------------------------
883
884   void CDomain::checkMask(void)
885   {
886      if (!mask_1d.isEmpty() && !mask_2d.isEmpty())
887        ERROR("CDomain::checkMask(void)",
888              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
889              << "Both mask_1d and mask_2d are defined but only one can be used at the same time." << std::endl
890              << "Please define only one mask: 'mask_1d' or 'mask_2d'.");
891
892      if (!mask_1d.isEmpty() && mask_2d.isEmpty())
893      {
894        if (mask_1d.numElements() != i_index.numElements())
895          ERROR("CDomain::checkMask(void)",
896                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
897                << "'mask_1d' does not have the same size as the local domain." << std::endl
898                << "Local size is " << i_index.numElements() << "." << std::endl
899                << "Mask size is " << mask_1d.numElements() << ".");
900      }
901
902      if (mask_1d.isEmpty() && !mask_2d.isEmpty())
903      {
904        if (mask_2d.extent(0) != ni || mask_2d.extent(1) != nj)
905          ERROR("CDomain::checkMask(void)",
906                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
907                << "The mask does not have the same size as the local domain." << std::endl
908                << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
909                << "Mask size is " << mask_2d.extent(0) << " x " << mask_2d.extent(1) << ".");
910      }
911
912      if (!mask_2d.isEmpty())
913      {
914        mask_1d.resize(mask_2d.extent(0) * mask_2d.extent(1));
915        for (int j = 0; j < nj; ++j)
916          for (int i = 0; i < ni; ++i) mask_1d(i+j*ni) = mask_2d(i,j);
917        mask_2d.reset();
918      }
919      else if (mask_1d.isEmpty())
920      {
921        mask_1d.resize(i_index.numElements());
922        for (int i = 0; i < i_index.numElements(); ++i) mask_1d(i) = true;
923      }
924   }
925
926   //----------------------------------------------------------------
927
928   void CDomain::checkDomainData(void)
929   {
930      if (data_dim.isEmpty())
931      {
932        data_dim.setValue(1);
933      }
934      else if (!(data_dim.getValue() == 1 || data_dim.getValue() == 2))
935      {
936        ERROR("CDomain::checkDomainData(void)",
937              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
938              << "The data dimension is invalid, 'data_dim' must be 1 or 2 not << " << data_dim.getValue() << ".");
939      }
940
941      if (data_ibegin.isEmpty())
942         data_ibegin.setValue(0);
943      if (data_jbegin.isEmpty())
944         data_jbegin.setValue(0);
945
946      if (data_ni.isEmpty())
947      {
948        data_ni.setValue((data_dim == 1) ? (ni.getValue() * nj.getValue()) : ni.getValue());
949      }
950      else if (data_ni.getValue() < 0)
951      {
952        ERROR("CDomain::checkDomainData(void)",
953              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
954              << "The data size cannot be negative ('data_ni' = " << data_ni.getValue() << ").");
955      }
956
957      if (data_nj.isEmpty())
958      {
959        data_nj.setValue((data_dim.getValue() == 1) ? (ni.getValue() * nj.getValue()) : nj.getValue());
960      }
961      else if (data_nj.getValue() < 0)
962      {
963        ERROR("CDomain::checkDomainData(void)",
964              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
965              << "The data size cannot be negative ('data_nj' = " << data_nj.getValue() << ").");
966      }
967   }
968
969   //----------------------------------------------------------------
970
971   void CDomain::checkCompression(void)
972   {
973      if (!data_i_index.isEmpty())
974      {
975        if (!data_j_index.isEmpty() &&
976            data_j_index.numElements() != data_i_index.numElements())
977        {
978           ERROR("CDomain::checkCompression(void)",
979                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
980                 << "'data_i_index' and 'data_j_index' arrays must have the same size." << std::endl
981                 << "'data_i_index' size = " << data_i_index.numElements() << std::endl
982                 << "'data_j_index' size = " << data_j_index.numElements());
983        }
984
985        if (2 == data_dim)
986        {
987          if (data_j_index.isEmpty())
988          {
989             ERROR("CDomain::checkCompression(void)",
990                   << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
991                   << "'data_j_index' must be defined when 'data_i_index' is set and 'data_dim' is 2.");
992          }
993        }
994        else // (1 == data_dim)
995        {
996          if (data_j_index.isEmpty())
997          {
998            data_j_index.resize(data_ni);
999            for (int j = 0; j < data_ni; ++j) data_j_index(j) = 0;
1000          }
1001        }
1002      }
1003      else
1004      {
1005        if (data_dim == 2 && !data_j_index.isEmpty())
1006          ERROR("CDomain::checkCompression(void)",
1007                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1008                << "'data_i_index' must be defined when 'data_j_index' is set and 'data_dim' is 2.");
1009
1010        if (1 == data_dim)
1011        {
1012          data_i_index.resize(data_ni);
1013          data_j_index.resize(data_ni);
1014
1015          for (int i = 0; i < data_ni; ++i)
1016          {
1017            data_i_index(i) = i;
1018            data_j_index(i) = 0;
1019          }
1020        }
1021        else // (data_dim == 2)
1022        {
1023          const int dsize = data_ni * data_nj;
1024          data_i_index.resize(dsize);
1025          data_j_index.resize(dsize);
1026
1027          for(int count = 0, j = 0; j < data_nj; ++j)
1028          {
1029            for(int i = 0; i < data_ni; ++i, ++count)
1030            {
1031              data_i_index(count) = i;
1032              data_j_index(count) = j;
1033            }
1034          }
1035        }
1036      }
1037   }
1038
1039   //----------------------------------------------------------------
1040   void CDomain::computeLocalMask(void)
1041   {
1042     localMask.resize(ni*nj) ;
1043     localMask=false ;
1044     size_t zoom_ibegin= global_zoom_ibegin ;
1045     size_t zoom_iend= global_zoom_ibegin+global_zoom_ni-1 ;
1046     size_t zoom_jbegin= global_zoom_jbegin ;
1047     size_t zoom_jend= global_zoom_jbegin+global_zoom_nj-1 ;
1048
1049
1050     size_t dn=data_i_index.numElements() ;
1051     int i,j ;
1052     size_t k,ind ;
1053
1054     for(k=0;k<dn;k++)
1055     {
1056       if (data_dim==2)
1057       {
1058          i=data_i_index(k)+data_ibegin ;
1059          j=data_j_index(k)+data_jbegin ;
1060       }
1061       else
1062       {
1063          i=(data_i_index(k)+data_ibegin)%ni ;
1064          j=(data_i_index(k)+data_ibegin)/ni ;
1065       }
1066
1067       if (i>=0 && i<ni && j>=0 && j<nj)
1068         if (i+ibegin>=zoom_ibegin && i+ibegin<=zoom_iend && j+jbegin>=zoom_jbegin && j+jbegin<=zoom_jend)
1069         {
1070           ind=i+ni*j ;
1071           localMask(ind)=mask_1d(ind) ;
1072         }
1073     }
1074   }
1075
1076   void CDomain::checkEligibilityForCompressedOutput(void)
1077   {
1078     // We don't check if the mask or the indexes are valid here, just if they have been defined at this point.
1079     isCompressible_ = !mask_1d.isEmpty() || !mask_2d.isEmpty() || !data_i_index.isEmpty();
1080   }
1081
1082   //----------------------------------------------------------------
1083
1084   /*
1085     Fill in longitude and latitude value from clients (or models) into internal values lonvalue, latvalue which
1086     will be used by XIOS.
1087   */
1088   void CDomain::completeLonLatClient(void)
1089   {
1090     bool lonlatValueExisted = (0 != lonvalue.numElements()) || (0 != latvalue.numElements());
1091     if (!lonvalue_2d.isEmpty() && !lonlatValueExisted)
1092     {
1093       lonvalue.resize(ni * nj);
1094       latvalue.resize(ni * nj);
1095       if (hasBounds)
1096       {
1097         bounds_lonvalue.resize(nvertex, ni * nj);
1098         bounds_latvalue.resize(nvertex, ni * nj);
1099       }
1100
1101       for (int j = 0; j < nj; ++j)
1102       {
1103         for (int i = 0; i < ni; ++i)
1104         {
1105           int k = j * ni + i;
1106
1107           lonvalue(k) = lonvalue_2d(i,j);
1108           latvalue(k) = latvalue_2d(i,j);
1109
1110           if (hasBounds)
1111           {
1112             for (int n = 0; n < nvertex; ++n)
1113             {
1114               bounds_lonvalue(n,k) = bounds_lon_2d(n,i,j);
1115               bounds_latvalue(n,k) = bounds_lat_2d(n,i,j);
1116             }
1117           }
1118         }
1119       }
1120     }
1121     else if (!lonvalue_1d.isEmpty()  && !lonlatValueExisted)
1122     {
1123       if (type_attr::rectilinear == type)
1124       {
1125         if (ni == lonvalue_1d.numElements() && nj == latvalue_1d.numElements())
1126         {
1127           lonvalue.resize(ni * nj);
1128           latvalue.resize(ni * nj);
1129           if (hasBounds)
1130           {
1131             bounds_lonvalue.resize(nvertex, ni * nj);
1132             bounds_latvalue.resize(nvertex, ni * nj);
1133           }
1134
1135           for (int j = 0; j < nj; ++j)
1136           {
1137             for (int i = 0; i < ni; ++i)
1138             {
1139               int k = j * ni + i;
1140
1141               lonvalue(k) = lonvalue_1d(i);
1142               latvalue(k) = latvalue_1d(j);
1143
1144               if (hasBounds)
1145               {
1146                 for (int n = 0; n < nvertex; ++n)
1147                 {
1148                   bounds_lonvalue(n,k) = bounds_lon_1d(n,i);
1149                   bounds_latvalue(n,k) = bounds_lat_1d(n,j);
1150                 }
1151               }
1152             }
1153           }
1154         }
1155         else if (i_index.numElements() == lonvalue_1d.numElements() && j_index.numElements() == latvalue_1d.numElements()  && !lonlatValueExisted)
1156         {
1157           lonvalue.reference(lonvalue_1d);
1158           latvalue.reference(latvalue_1d);
1159            if (hasBounds)
1160           {
1161             bounds_lonvalue.reference(bounds_lon_1d);
1162             bounds_latvalue.reference(bounds_lat_1d);
1163           }
1164         }
1165         else
1166           ERROR("CDomain::completeLonClient(void)",
1167                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1168                 << "'lonvalue_1d' and 'latvalue_1d' does not have the same size as the local domain." << std::endl
1169                 << "'lonvalue_1d' size is " << lonvalue_1d.numElements() 
1170                 << " and 'latvalue_1d' size is " << latvalue_1d.numElements() << std::endl
1171                 << " They should be correspondingly " << ni.getValue() << " and "  << nj.getValue() << " or " << std::endl
1172                 << i_index.numElements() << " and "  << j_index.numElements() << ".");
1173       }
1174       else if (type == type_attr::curvilinear || type == type_attr::unstructured  && !lonlatValueExisted)
1175       {
1176         lonvalue.reference(lonvalue_1d);
1177         latvalue.reference(latvalue_1d);
1178         if (hasBounds)
1179         {
1180           bounds_lonvalue.reference(bounds_lon_1d);
1181           bounds_latvalue.reference(bounds_lat_1d);
1182         }
1183       }
1184     }
1185   }
1186
1187   /*
1188     Convert internal longitude latitude value used by XIOS to "lonvalue_*" which can be retrieved with Fortran interface
1189   */
1190   void CDomain::convertLonLatValue(void)
1191   {
1192     bool lonlatValueExisted = (0 != lonvalue.numElements()) || (0 != latvalue.numElements());
1193     if (!lonvalue_2d.isEmpty() && lonlatValueExisted)
1194     {
1195       lonvalue_2d.resize(ni,nj);
1196       latvalue_2d.resize(ni,nj);
1197       if (hasBounds)
1198       {
1199         bounds_lon_2d.resize(nvertex, ni, nj);
1200         bounds_lat_2d.resize(nvertex, ni, nj);
1201       }
1202
1203       for (int j = 0; j < nj; ++j)
1204       {
1205         for (int i = 0; i < ni; ++i)
1206         {
1207           int k = j * ni + i;
1208
1209           lonvalue_2d(i,j) = lonvalue(k);
1210           latvalue_2d(i,j) = latvalue(k);
1211
1212           if (hasBounds)
1213           {
1214             for (int n = 0; n < nvertex; ++n)
1215             {
1216               bounds_lon_2d(n,i,j) = bounds_lonvalue(n,k);
1217               bounds_lat_2d(n,i,j) = bounds_latvalue(n,k);
1218             }
1219           }
1220         }
1221       }
1222     }
1223     else if (!lonvalue_1d.isEmpty()  && lonlatValueExisted)
1224     {
1225       if (type_attr::rectilinear == type)
1226       {
1227         if (ni == lonvalue_1d.numElements() && nj == latvalue_1d.numElements())
1228         {
1229           lonvalue.resize(ni * nj);
1230           latvalue.resize(ni * nj);
1231           if (hasBounds)
1232           {
1233             bounds_lonvalue.resize(nvertex, ni * nj);
1234             bounds_latvalue.resize(nvertex, ni * nj);
1235           }
1236
1237           for (int j = 0; j < nj; ++j)
1238           {
1239             for (int i = 0; i < ni; ++i)
1240             {
1241               int k = j * ni + i;
1242
1243               lonvalue(k) = lonvalue_1d(i);
1244               latvalue(k) = latvalue_1d(j);
1245
1246               if (hasBounds)
1247               {
1248                 for (int n = 0; n < nvertex; ++n)
1249                 {
1250                   bounds_lonvalue(n,k) = bounds_lon_1d(n,i);
1251                   bounds_latvalue(n,k) = bounds_lat_1d(n,j);
1252                 }
1253               }
1254             }
1255           }
1256         }
1257         else if (i_index.numElements() == lonvalue_1d.numElements() && j_index.numElements() == latvalue_1d.numElements()  && !lonlatValueExisted)
1258         {
1259           lonvalue.reference(lonvalue_1d);
1260           latvalue.reference(latvalue_1d);
1261            if (hasBounds)
1262           {
1263             bounds_lonvalue.reference(bounds_lon_1d);
1264             bounds_latvalue.reference(bounds_lat_1d);
1265           }
1266         }
1267         else
1268           ERROR("CDomain::completeLonClient(void)",
1269                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1270                 << "'lonvalue_1d' and 'latvalue_1d' does not have the same size as the local domain." << std::endl
1271                 << "'lonvalue_1d' size is " << lonvalue_1d.numElements() 
1272                 << " and 'latvalue_1d' size is " << latvalue_1d.numElements() << std::endl
1273                 << " They should be correspondingly " << ni.getValue() << " and "  << nj.getValue() << " or " << std::endl
1274                 << i_index.numElements() << " and "  << j_index.numElements() << ".");
1275       }
1276       else if (type == type_attr::curvilinear || type == type_attr::unstructured  && !lonlatValueExisted)
1277       {
1278         lonvalue.reference(lonvalue_1d);
1279         latvalue.reference(latvalue_1d);
1280         if (hasBounds)
1281         {
1282           bounds_lonvalue.reference(bounds_lon_1d);
1283           bounds_latvalue.reference(bounds_lat_1d);
1284         }
1285       }
1286     }
1287   }
1288
1289
1290   void CDomain::checkBounds(void)
1291   {
1292     bool hasBoundValues = (0 != bounds_lonvalue.numElements()) || (0 != bounds_latvalue.numElements());
1293     if (!nvertex.isEmpty() && nvertex > 0 && !hasBoundValues)
1294     {
1295       if (!bounds_lon_1d.isEmpty() && !bounds_lon_2d.isEmpty())
1296         ERROR("CDomain::checkBounds(void)",
1297               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1298               << "Only one longitude boundary attribute can be used but both 'bounds_lon_1d' and 'bounds_lon_2d' are defined." << std::endl
1299               << "Define only one longitude boundary attribute: 'bounds_lon_1d' or 'bounds_lon_2d'.");
1300
1301       if (!bounds_lat_1d.isEmpty() && !bounds_lat_2d.isEmpty())
1302         ERROR("CDomain::checkBounds(void)",
1303               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1304               << "Only one latitude boundary attribute can be used but both 'bounds_lat_1d' and 'bounds_lat_2d' are defined." << std::endl
1305               << "Define only one latitude boundary attribute: 'bounds_lat_1d' or 'bounds_lat_2d'.");
1306
1307       if ((!bounds_lon_1d.isEmpty() && bounds_lat_1d.isEmpty()) || (bounds_lon_1d.isEmpty() && !bounds_lat_1d.isEmpty()))
1308       {
1309         ERROR("CDomain::checkBounds(void)",
1310               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1311               << "Only 'bounds_lon_1d' or 'bounds_lat_1d' is defined." << std::endl
1312               << "Please define either both attributes or none.");
1313       }
1314
1315       if ((!bounds_lon_2d.isEmpty() && bounds_lat_2d.isEmpty()) || (bounds_lon_2d.isEmpty() && !bounds_lat_2d.isEmpty()))
1316       {
1317         ERROR("CDomain::checkBounds(void)",
1318               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1319               << "Only 'bounds_lon_2d' or 'bounds_lat_2d' is defined." << std::endl
1320               << "Please define either both attributes or none.");
1321       }
1322
1323       if (!bounds_lon_1d.isEmpty() && nvertex.getValue() != bounds_lon_1d.extent(0))
1324         ERROR("CDomain::checkBounds(void)",
1325               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1326               << "'bounds_lon_1d' dimension is not compatible with 'nvertex'." << std::endl
1327               << "'bounds_lon_1d' dimension is " << bounds_lon_1d.extent(1)
1328               << " but nvertex is " << nvertex.getValue() << ".");
1329
1330       if (!bounds_lon_2d.isEmpty() && nvertex.getValue() != bounds_lon_2d.extent(0))
1331         ERROR("CDomain::checkBounds(void)",
1332               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1333               << "'bounds_lon_2d' dimension is not compatible with 'nvertex'." << std::endl
1334               << "'bounds_lon_2d' dimension is " << bounds_lon_2d.extent(2)
1335               << " but nvertex is " << nvertex.getValue() << ".");
1336
1337       if (!bounds_lon_1d.isEmpty() && lonvalue_1d.isEmpty())
1338         ERROR("CDomain::checkBounds(void)",
1339               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1340               << "Since 'bounds_lon_1d' is defined, 'lonvalue_1d' must be defined too." << std::endl);
1341
1342       if (!bounds_lon_2d.isEmpty() && lonvalue_2d.isEmpty())
1343         ERROR("CDomain::checkBounds(void)",
1344               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1345               << "Since 'bounds_lon_2d' is defined, 'lonvalue_2d' must be defined too." << std::endl);
1346
1347       if (!bounds_lat_1d.isEmpty() && nvertex.getValue() != bounds_lat_1d.extent(0))
1348         ERROR("CDomain::checkBounds(void)",
1349               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1350               << "'bounds_lat_1d' dimension is not compatible with 'nvertex'." << std::endl
1351               << "'bounds_lat_1d' dimension is " << bounds_lat_1d.extent(1)
1352               << " but nvertex is " << nvertex.getValue() << ".");
1353
1354       if (!bounds_lat_2d.isEmpty() && nvertex.getValue() != bounds_lat_2d.extent(0))
1355         ERROR("CDomain::checkBounds(void)",
1356               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1357               << "'bounds_lat_2d' dimension is not compatible with 'nvertex'." << std::endl
1358               << "'bounds_lat_2d' dimension is " << bounds_lat_2d.extent(2)
1359               << " but nvertex is " << nvertex.getValue() << ".");
1360
1361       if (!bounds_lat_1d.isEmpty() && latvalue_1d.isEmpty())
1362         ERROR("CDomain::checkBounds(void)",
1363               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1364               << "Since 'bounds_lat_1d' is defined, 'latvalue_1d' must be defined too." << std::endl);
1365
1366       if (!bounds_lat_2d.isEmpty() && latvalue_2d.isEmpty())
1367         ERROR("CDomain::checkBounds(void)",
1368               << "Since 'bounds_lat_2d' is defined, 'latvalue_2d' must be defined too." << std::endl);
1369
1370       hasBounds = true;
1371     }
1372     else if (hasBoundValues)
1373     {
1374       hasBounds = true;       
1375     }
1376     else
1377     {
1378       hasBounds = false;
1379       nvertex = 0;
1380     }
1381   }
1382
1383   void CDomain::checkArea(void)
1384   {
1385     bool hasAreaValue = (0 != areavalue.numElements());
1386     hasArea = !area.isEmpty() || !areavalue.isEmpty();
1387     if (hasArea)
1388     {
1389       if (area.extent(0) != ni || area.extent(1) != nj)
1390       {
1391         ERROR("CDomain::checkArea(void)",
1392               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1393               << "The area does not have the same size as the local domain." << std::endl
1394               << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
1395               << "Area size is " << area.extent(0) << " x " << area.extent(1) << ".");
1396       }
1397       if (areavalue.isEmpty())
1398       {
1399          areavalue.resize(ni*nj);
1400         for (int j = 0; j < nj; ++j)
1401         {
1402           for (int i = 0; i < ni; ++i)
1403           {
1404             int k = j * ni + i;
1405             areavalue(k) = area(i,j);
1406           }
1407         }
1408       }
1409     }
1410   }
1411
1412   void CDomain::checkLonLat()
1413   {
1414     hasLonLat = (!latvalue_1d.isEmpty() && !lonvalue_1d.isEmpty()) ||
1415                 (!latvalue_2d.isEmpty() && !lonvalue_2d.isEmpty());
1416     bool hasLonLatValue = (0 != lonvalue.numElements()) || (0 != latvalue.numElements());
1417     if (hasLonLat && !hasLonLatValue)
1418     {
1419       if (!lonvalue_1d.isEmpty() && !lonvalue_2d.isEmpty())
1420         ERROR("CDomain::checkLonLat()",
1421               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1422               << "Only one longitude attribute can be used but both 'lonvalue_1d' and 'lonvalue_2d' are defined." << std::endl
1423               << "Define only one longitude attribute: 'lonvalue_1d' or 'lonvalue_2d'.");
1424
1425       if (!lonvalue_1d.isEmpty() && lonvalue_2d.isEmpty())
1426       {
1427         if ((type_attr::rectilinear != type) && (lonvalue_1d.numElements() != i_index.numElements()))
1428           ERROR("CDomain::checkLonLat()",
1429                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1430                 << "'lonvalue_1d' does not have the same size as the local domain." << std::endl
1431                 << "Local size is " << i_index.numElements() << "." << std::endl
1432                 << "'lonvalue_1d' size is " << lonvalue_1d.numElements() << ".");
1433       }
1434
1435       if (lonvalue_1d.isEmpty() && !lonvalue_2d.isEmpty())
1436       {
1437         if (lonvalue_2d.extent(0) != ni || lonvalue_2d.extent(1) != nj)
1438           ERROR("CDomain::checkLonLat()",
1439                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1440                 << "'lonvalue_2d' does not have the same size as the local domain." << std::endl
1441                 << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
1442                 << "'lonvalue_2d' size is " << lonvalue_2d.extent(0) << " x " << lonvalue_2d.extent(1) << ".");
1443       }
1444
1445       if (!latvalue_1d.isEmpty() && !latvalue_2d.isEmpty())
1446         ERROR("CDomain::checkLonLat()",
1447               << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1448               << "Only one latitude attribute can be used but both 'latvalue_1d' and 'latvalue_2d' are defined." << std::endl
1449               << "Define only one latitude attribute: 'latvalue_1d' or 'latvalue_2d'.");
1450
1451       if (!latvalue_1d.isEmpty() && latvalue_2d.isEmpty())
1452       {
1453         if ((type_attr::rectilinear != type) && (latvalue_1d.numElements() != i_index.numElements()))
1454           ERROR("CDomain::checkLonLat()",
1455                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1456                 << "'latvalue_1d' does not have the same size as the local domain." << std::endl
1457                 << "Local size is " << i_index.numElements() << "." << std::endl
1458                 << "'latvalue_1d' size is " << latvalue_1d.numElements() << ".");
1459       }
1460
1461       if (latvalue_1d.isEmpty() && !latvalue_2d.isEmpty())
1462       {
1463         if (latvalue_2d.extent(0) != ni || latvalue_2d.extent(1) != nj)
1464           ERROR("CDomain::checkLonLat()",
1465                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
1466                 << "'latvalue_2d' does not have the same size as the local domain." << std::endl
1467                 << "Local size is " << ni.getValue() << " x " << nj.getValue() << "." << std::endl
1468                 << "'latvalue_2d' size is " << latvalue_2d.extent(0) << " x " << latvalue_2d.extent(1) << ".");
1469       }
1470     }
1471   }
1472
1473   void CDomain::checkAttributesOnClientAfterTransformation()
1474   {
1475     CContext* context=CContext::getCurrent() ;
1476
1477     if (this->isClientAfterTransformationChecked) return;
1478     if (context->hasClient)
1479     {
1480       // this->checkMask();
1481      this->computeConnectedClients();
1482       // if (hasLonLat || hasArea || isCompressible_) this->computeConnectedClients();
1483       if (hasLonLat)
1484         if (!context->hasServer)
1485           this->completeLonLatClient();
1486     }
1487
1488     this->isClientAfterTransformationChecked = true;
1489   }
1490
1491   //----------------------------------------------------------------
1492   // Divide function checkAttributes into 2 seperate ones
1493   // This function only checks all attributes of current domain
1494   void CDomain::checkAttributesOnClient()
1495   {
1496     if (this->isClientChecked) return;
1497     CContext* context=CContext::getCurrent();
1498
1499      if (context->hasClient && !context->hasServer)
1500      {
1501        this->checkDomain();
1502        this->checkBounds();
1503        this->checkArea();
1504        this->checkLonLat();
1505      }
1506
1507      if (context->hasClient && !context->hasServer)
1508      { // Ct client uniquement
1509         this->checkMask();
1510         this->checkDomainData();
1511         this->checkCompression();
1512         this->computeLocalMask() ;
1513      }
1514      else
1515      { // Ct serveur uniquement
1516      }
1517
1518      this->isClientChecked = true;
1519   }
1520
1521   // Send all checked attributes to server
1522   void CDomain::sendCheckedAttributes()
1523   {
1524     if (!this->isClientChecked) checkAttributesOnClient();
1525     if (!this->isClientAfterTransformationChecked) checkAttributesOnClientAfterTransformation();
1526     CContext* context=CContext::getCurrent() ;
1527
1528     if (this->isChecked) return;
1529     if (context->hasClient)
1530     {
1531//       this->completeLonLatClient();
1532       sendAttributes();
1533     }
1534     this->isChecked = true;
1535   }
1536
1537   void CDomain::checkAttributes(void)
1538   {
1539      if (this->isChecked) return;
1540      CContext* context=CContext::getCurrent() ;
1541
1542      this->checkDomain();
1543      this->checkLonLat();
1544      this->checkBounds();
1545      this->checkArea();
1546
1547      if (context->hasClient)
1548      { // Ct client uniquement
1549         this->checkMask();
1550         this->checkDomainData();
1551         this->checkCompression();
1552         this->computeLocalMask() ;
1553
1554      }
1555      else
1556      { // Ct serveur uniquement
1557      }
1558
1559      if (context->hasClient)
1560      {
1561        this->computeConnectedClients();
1562        this->completeLonLatClient();
1563      }
1564
1565      this->isChecked = true;
1566   }
1567
1568
1569  // void CDomain::computeConnectedClients(const std::vector<int>& globalDim, int orderPositionInGrid,
1570  //                                     CServerDistributionDescription::ServerDistributionType distType)
1571  /*!
1572     Compute the connection of a client to other clients to determine which clients to send attributes to
1573  */
1574  void CDomain::computeConnectedClients()
1575  {
1576    CContext* context=CContext::getCurrent() ;
1577
1578    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
1579    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
1580    for (int p = 0; p < nbSrvPools; ++p)
1581    {
1582      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
1583      int nbServer = client->serverSize;
1584      int rank     = client->clientRank;
1585      bool doComputeGlobalIndexServer = true;
1586
1587      int i,j,i_ind,j_ind, nbIndex, nbIndexZoom;
1588      int global_zoom_iend=global_zoom_ibegin+global_zoom_ni-1;
1589      int global_zoom_jend=global_zoom_jbegin+global_zoom_nj-1;
1590
1591      // Precompute number of index
1592      int globalIndexCountZoom = 0;
1593      nbIndex = i_index.numElements();
1594
1595      if (doZoomByIndex_) 
1596      {
1597        globalIndexCountZoom = zoom_i_index.numElements();
1598      }
1599      else 
1600      {
1601        for (i = 0; i < nbIndex; ++i)
1602        {
1603          i_ind=i_index(i);
1604          j_ind=j_index(i);
1605
1606          if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1607          {
1608            ++globalIndexCountZoom;
1609          }
1610        }
1611      }
1612
1613      // for (i = 0; i < nbIndex; ++i)
1614      // {
1615      //   i_ind=i_index(i);
1616      //   j_ind=j_index(i);
1617
1618      //   if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1619      //   {
1620      //     ++globalIndexCountZoom;
1621      //   }
1622      // }
1623
1624      // int globalIndexWrittenCount = 0;
1625      // if (isCompressible_)
1626      // {
1627      //   for (i = 0; i < data_i_index.numElements(); ++i)
1628      //   {
1629      //     i_ind = CDistributionClient::getDomainIndex(data_i_index(i), data_j_index(i),
1630      //                                                 data_ibegin, data_jbegin, data_dim, ni,
1631      //                                                 j_ind);
1632      //     if (i_ind >= 0 && i_ind < ni && j_ind >= 0 && j_ind < nj && mask_1d(i_ind + j_ind * ni))
1633      //     {
1634      //       i_ind += ibegin;
1635      //       j_ind += jbegin;
1636      //       if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1637      //         ++globalIndexWrittenCount;
1638      //     }
1639      //   }
1640      // }
1641
1642      // Fill in index
1643
1644      CArray<size_t,1> globalIndexDomainZoom(globalIndexCountZoom);
1645      CArray<size_t,1> localIndexDomainZoom(globalIndexCountZoom);
1646      CArray<size_t,1> globalIndexDomain(nbIndex);
1647      size_t globalIndex;
1648      int globalIndexCount = 0;
1649
1650      for (i = 0; i < nbIndex; ++i)
1651      {
1652        i_ind=i_index(i);
1653        j_ind=j_index(i);
1654        globalIndex = i_ind + j_ind * ni_glo;
1655        globalIndexDomain(i) = globalIndex;               
1656      }
1657
1658      if (globalLocalIndexMap_.empty())
1659      {
1660        for (i = 0; i < nbIndex; ++i)
1661          globalLocalIndexMap_[globalIndexDomain(i)] = i;
1662      }
1663
1664      globalIndexCountZoom = 0;
1665      if (doZoomByIndex_) 
1666      {
1667        int nbIndexZoom = zoom_i_index.numElements();       
1668       
1669        for (i = 0; i < nbIndexZoom; ++i)
1670        {
1671          i_ind=zoom_i_index(i);
1672          j_ind=zoom_j_index(i);
1673          globalIndex = i_ind + j_ind * ni_glo;
1674          globalIndexDomainZoom(globalIndexCountZoom) = globalIndex;
1675          ++globalIndexCountZoom;
1676          // if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1677          // {
1678          //   globalIndexDomainZoom(globalIndexCountZoom) = globalIndex;
1679          //   localIndexDomainZoom(globalIndexCountZoom) = i;
1680          //   ++globalIndexCountZoom;
1681          // }
1682        }
1683      }
1684      else 
1685      {
1686          int global_zoom_iend=global_zoom_ibegin+global_zoom_ni-1;
1687          int global_zoom_jend=global_zoom_jbegin+global_zoom_nj-1;
1688          for (i = 0; i < nbIndex; ++i)
1689          {
1690            i_ind=i_index(i);
1691            j_ind=j_index(i);
1692            globalIndex = i_ind + j_ind * ni_glo;
1693            if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1694            {
1695              globalIndexDomainZoom(globalIndexCountZoom) = globalIndex;
1696              ++globalIndexCountZoom;
1697            }
1698          }
1699
1700          int iend = ibegin + ni -1;
1701          int jend = jbegin + nj -1;
1702          zoom_ibegin = global_zoom_ibegin > ibegin ? global_zoom_ibegin : ibegin;
1703          int zoom_iend   = global_zoom_iend < iend ? zoom_iend : iend ;
1704          zoom_ni     = zoom_iend-zoom_ibegin+1 ;
1705
1706          zoom_jbegin = global_zoom_jbegin > jbegin ? global_zoom_jbegin : jbegin ;
1707          int zoom_jend   = global_zoom_jend < jend ? zoom_jend : jend;
1708          zoom_nj     = zoom_jend-zoom_jbegin+1;
1709      }
1710
1711      // if (globalLocalIndexZoomMap_.empty())
1712      // {
1713      //   int nbIndexZoom = globalIndexDomainZoom.numElements();
1714      //   for (i = 0; i < nbIndex; ++i)
1715      //     globalLocalIndexZoomMap_[globalIndexDomainZoom(i)] = i;
1716      // }
1717
1718
1719      // CArray<int,1> globalIndexWrittenDomain(globalIndexWrittenCount);
1720      // if (isCompressible_)
1721      // {
1722      //   globalIndexWrittenCount = 0;
1723      //   for (i = 0; i < data_i_index.numElements(); ++i)
1724      //   {
1725      //     i_ind = CDistributionClient::getDomainIndex(data_i_index(i), data_j_index(i),
1726      //                                                 data_ibegin, data_jbegin, data_dim, ni,
1727      //                                                 j_ind);
1728      //     if (i_ind >= 0 && i_ind < ni && j_ind >= 0 && j_ind < nj && mask_1d(i_ind + j_ind * ni))
1729      //     {
1730      //       i_ind += ibegin;
1731      //       j_ind += jbegin;
1732      //       if (i_ind >= global_zoom_ibegin && i_ind <= global_zoom_iend && j_ind >= global_zoom_jbegin && j_ind <= global_zoom_jend)
1733      //       {
1734      //         globalIndexWrittenDomain(globalIndexWrittenCount) = i_ind + j_ind * ni_glo;
1735      //         ++globalIndexWrittenCount;
1736      //       }
1737      //     }
1738      //   }
1739      // }
1740
1741      size_t globalSizeIndex = 1, indexBegin, indexEnd;
1742      int range, clientSize = client->clientSize;
1743      std::vector<int> nGlobDomain(2);
1744      nGlobDomain[0] = this->ni_glo;
1745      nGlobDomain[1] = this->nj_glo;
1746      for (int i = 0; i < nGlobDomain.size(); ++i) globalSizeIndex *= nGlobDomain[i];
1747      indexBegin = 0;
1748      if (globalSizeIndex <= clientSize)
1749      {
1750        indexBegin = rank%globalSizeIndex;
1751        indexEnd = indexBegin;
1752      }
1753      else
1754      {
1755        for (int i = 0; i < clientSize; ++i)
1756        {
1757          range = globalSizeIndex / clientSize;
1758          if (i < (globalSizeIndex%clientSize)) ++range;
1759          if (i == client->clientRank) break;
1760          indexBegin += range;
1761        }
1762        indexEnd = indexBegin + range - 1;
1763      }
1764
1765      CServerDistributionDescription serverDescription(nGlobDomain, nbServer);
1766      if (isUnstructed_) serverDescription.computeServerGlobalIndexInRange(std::make_pair<size_t,size_t>(indexBegin, indexEnd), 0);
1767      else serverDescription.computeServerGlobalIndexInRange(std::make_pair<size_t,size_t>(indexBegin, indexEnd), 1);
1768
1769      CClientServerMapping* clientServerMap = new CClientServerMappingDistributed(serverDescription.getGlobalIndexRange(),
1770                                                                                  client->intraComm);
1771      clientServerMap->computeServerIndexMapping(globalIndexDomain);
1772      CClientServerMapping::GlobalIndexMap& globalIndexDomainOnServer = clientServerMap->getGlobalIndexOnServer();
1773
1774      CClientServerMapping::GlobalIndexMap::const_iterator it  = globalIndexDomainOnServer.begin(),
1775                                                           ite = globalIndexDomainOnServer.end();
1776      // typedef XIOSBinarySearchWithIndex<size_t> BinarySearch;
1777      // std::vector<int>::iterator itVec;
1778
1779      // indSrv_.clear();
1780      // indWrittenSrv_.clear();
1781      // for (; it != ite; ++it)
1782      // {
1783      //   int rank = it->first;
1784      //   int indexSize = it->second.size();
1785      //   std::vector<int> permutIndex(indexSize);
1786      //   XIOSAlgorithms::fillInIndex(indexSize, permutIndex);
1787      //   XIOSAlgorithms::sortWithIndex<size_t, CVectorStorage>(it->second, permutIndex);
1788      //   BinarySearch binSearch(it->second);
1789      //   int nb = globalIndexDomainZoom.numElements();
1790      //   for (int i = 0; i < nb; ++i)
1791      //   {
1792      //     if (binSearch.search(permutIndex.begin(), permutIndex.end(), globalIndexDomainZoom(i), itVec))
1793      //     {
1794      //       indSrv_[rank].push_back(localIndexDomainZoom(i));
1795      //     }
1796      //   }
1797      //   for (int i = 0; i < globalIndexWrittenDomain.numElements(); ++i)
1798      //   {
1799      //     if (binSearch.search(permutIndex.begin(), permutIndex.end(), globalIndexWrittenDomain(i), itVec))
1800      //     {
1801      //       indWrittenSrv_[rank].push_back(globalIndexWrittenDomain(i));
1802      //     }
1803      //   }
1804      // }
1805
1806      connectedServerRank_.clear();
1807      for (it = globalIndexDomainOnServer.begin(); it != ite; ++it) {
1808        connectedServerRank_.push_back(it->first);
1809      }
1810
1811      indSrv_.swap(globalIndexDomainOnServer);
1812      nbConnectedClients_ = clientServerMap->computeConnectedClients(client->serverSize, client->clientSize, client->intraComm, connectedServerRank_);
1813
1814      clientServerMap->computeServerIndexMapping(globalIndexDomainZoom);
1815      CClientServerMapping::GlobalIndexMap& globalIndexDomainZoomOnServer = clientServerMap->getGlobalIndexOnServer();
1816      indZoomSrv_.swap(globalIndexDomainZoomOnServer);
1817     
1818     for (it = indZoomSrv_.begin(); it != indZoomSrv_.end(); ++it)
1819       connectedServerZoomRank_.push_back(it->first);
1820
1821      nbConnectedClientsZoom_ = clientServerMap->computeConnectedClients(client->serverSize, client->clientSize, client->intraComm, connectedServerZoomRank_);
1822
1823      delete clientServerMap;
1824    }
1825  }
1826
1827  const boost::unordered_map<int, vector<size_t> >& CDomain::getIndexServer() const
1828  {
1829    return indSrv_;
1830  }
1831
1832   void CDomain::computeWrittenIndex()
1833   { 
1834      if (computedWrittenIndex_) return;
1835      computedWrittenIndex_ = true;
1836
1837      CContext* context=CContext::getCurrent() ;
1838      CContextClient* client = context->client; 
1839
1840      std::vector<int> nBegin(2), nSize(2), nBeginGlobal(2), nGlob(2);
1841      nBegin[0]       = zoom_ibegin;  nBegin[1] = zoom_jbegin;
1842      nSize[0]        = zoom_ni;      nSize[1]  = zoom_nj;
1843      nBeginGlobal[0] = 0; nBeginGlobal[1] = 0;
1844      nGlob[0]        = ni_glo;   nGlob[1] = nj_glo;
1845      CDistributionServer srvDist(client->clientSize, nBegin, nSize, nBeginGlobal, nGlob); 
1846      const CArray<size_t,1>& writtenGlobalIndex  = srvDist.getGlobalIndex();
1847
1848      size_t nbWritten = 0, indGlo;     
1849      boost::unordered_map<size_t,size_t>::const_iterator itb = globalLocalIndexMap_.begin(),
1850                                                          ite = globalLocalIndexMap_.end(), it;         
1851      CArray<size_t,1>::const_iterator itSrvb = writtenGlobalIndex.begin(),
1852                                       itSrve = writtenGlobalIndex.end(), itSrv;
1853      localIndexToWriteOnServer.resize(writtenGlobalIndex.numElements());
1854      for (itSrv = itSrvb; itSrv != itSrve; ++itSrv)
1855      {
1856        indGlo = *itSrv;
1857        if (ite != globalLocalIndexMap_.find(indGlo))
1858        {
1859          localIndexToWriteOnServer(nbWritten) = globalLocalIndexMap_[indGlo];
1860          ++nbWritten;
1861        }                 
1862      }
1863   }
1864
1865  /*!
1866    Send all attributes from client to connected clients
1867    The attributes will be rebuilt on receiving side
1868  */
1869  void CDomain::sendAttributes()
1870  {
1871    sendDistributionAttributes();
1872    sendIndex();   
1873    // sendIndexZoom();
1874    sendMask();
1875    sendLonLat();
1876    sendArea();   
1877    sendDataIndex();
1878  }
1879
1880  /*!
1881    Send global index and zoom index from client to connected client(s)
1882    zoom index can be smaller than global index
1883  */
1884  void CDomain::sendIndex()
1885  {
1886    int ns, n, i, j, ind, nv, idx;
1887    CContext* context = CContext::getCurrent();
1888
1889    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
1890    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
1891    for (int p = 0; p < nbSrvPools; ++p)
1892    {
1893      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
1894
1895      CEventClient eventIndex(getType(), EVENT_ID_INDEX);
1896
1897      list<CMessage> list_msgsIndex;
1898      list<CArray<int,1> > list_indZoom, list_writtenInd, list_indGlob;
1899
1900      boost::unordered_map<int, vector<size_t> >::const_iterator itIndex, iteIndex;
1901      iteIndex = indSrv_.end();
1902      for (int k = 0; k < connectedServerRank_.size(); ++k)
1903      {
1904        int nbIndGlob = 0;
1905        int rank = connectedServerRank_[k];
1906        itIndex = indSrv_.find(rank);
1907        if (iteIndex != itIndex)
1908          nbIndGlob = itIndex->second.size();
1909
1910        list_indGlob.push_back(CArray<int,1>(nbIndGlob));       
1911
1912        CArray<int,1>& indGlob = list_indGlob.back();
1913        for (n = 0; n < nbIndGlob; ++n)
1914        {
1915          indGlob(n) = static_cast<int>(itIndex->second[n]);
1916        }
1917
1918        list_msgsIndex.push_back(CMessage());
1919        list_msgsIndex.back() << this->getId() << (int)type; // enum ne fonctionne pour les message => ToFix
1920        list_msgsIndex.back() << isCurvilinear;
1921        list_msgsIndex.back() << list_indGlob.back(); //list_indi.back() << list_indj.back();
1922       
1923        // if (isCompressible_)
1924        // {
1925        //   std::vector<int>& writtenIndSrc = indWrittenSrv_[rank];
1926        //   list_writtenInd.push_back(CArray<int,1>(writtenIndSrc.size()));
1927        //   CArray<int,1>& writtenInd = list_writtenInd.back();
1928
1929        //   for (n = 0; n < writtenInd.numElements(); ++n)
1930        //     writtenInd(n) = writtenIndSrc[n];
1931
1932        //   list_msgsIndex.back() << writtenInd;
1933        // }
1934
1935        eventIndex.push(rank, nbConnectedClients_[rank], list_msgsIndex.back());
1936      }
1937
1938      client->sendEvent(eventIndex);
1939    }
1940  }
1941
1942  /*!
1943    Send global index and zoom index from client to connected client(s)
1944    zoom index can be smaller than global index
1945  */
1946  void CDomain::sendIndexZoom()
1947  {
1948    int ns, n, i, j, ind, nv, idx;
1949    CContext* context = CContext::getCurrent();
1950
1951    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
1952    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
1953    for (int p = 0; p < nbSrvPools; ++p)
1954    {
1955      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
1956      CEventClient eventIndexZoom(getType(), EVENT_ID_INDEX_ZOOM);
1957
1958      list<CMessage> list_msgsIndex;
1959      list<CArray<int,1> > list_indZoom;
1960
1961      boost::unordered_map<int, vector<size_t> >::const_iterator itZoom, iteZoom;
1962      iteZoom = indZoomSrv_.end();
1963      for (int k = 0; k < connectedServerZoomRank_.size(); ++k)
1964      {
1965        int nbIndGlob = 0;
1966        int rank = connectedServerZoomRank_[k];
1967        int nbIndZoom = 0;
1968        itZoom = indZoomSrv_.find(rank);
1969        if (iteZoom != itZoom)
1970          nbIndZoom = itZoom->second.size();
1971       
1972        list_indZoom.push_back(CArray<int,1>(nbIndZoom));
1973        CArray<int,1>& indZoom = list_indZoom.back();
1974        for (n = 0; n < nbIndZoom; ++n)
1975        {
1976          indZoom(n) = static_cast<int>(itZoom->second[n]);
1977        }
1978
1979        list_msgsIndex.push_back(CMessage());
1980        list_msgsIndex.back() << this->getId(); // enum ne fonctionne pour les message => ToFix       
1981        list_msgsIndex.back() << list_indZoom.back() << doZoomByIndex_; //list_indi.back() << list_indj.back     
1982
1983        eventIndexZoom.push(rank, nbConnectedClientsZoom_[rank], list_msgsIndex.back());
1984      }
1985
1986      client->sendEvent(eventIndexZoom);
1987    }
1988  }
1989
1990  /*!
1991    Send distribution from client to other clients
1992    Because a client in a level knows correctly the grid distribution of client on the next level
1993    it calculates this distribution then sends it to the corresponding clients on the next level
1994  */
1995  void CDomain::sendDistributionAttributes(void)
1996  {
1997    CContext* context = CContext::getCurrent();
1998     // Use correct context client to send message
1999    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
2000    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
2001    for (int i = 0; i < nbSrvPools; ++i)
2002    {
2003      CContextClient* contextClientTmp = (context->hasServer) ? context->clientPrimServer[i]
2004                                                                         : context->client;   
2005      int nbServer = contextClientTmp->serverSize;
2006      std::vector<int> nGlobDomain(2);
2007      nGlobDomain[0] = this->ni_glo;
2008      nGlobDomain[1] = this->nj_glo;
2009
2010      CServerDistributionDescription serverDescription(nGlobDomain, nbServer);
2011      if (isUnstructed_) serverDescription.computeServerDistribution(false, 0);
2012      else serverDescription.computeServerDistribution(false, 1);
2013
2014      std::vector<std::vector<int> > serverIndexBegin = serverDescription.getServerIndexBegin();
2015      std::vector<std::vector<int> > serverDimensionSizes = serverDescription.getServerDimensionSizes();
2016
2017      CEventClient event(getType(),EVENT_ID_SERVER_ATTRIBUT);
2018      if (contextClientTmp->isServerLeader())
2019      {
2020        std::list<CMessage> msgs;
2021
2022        const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2023        for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2024        {
2025          // Use const int to ensure CMessage holds a copy of the value instead of just a reference
2026          const int ibegin_srv = serverIndexBegin[*itRank][0];
2027          const int jbegin_srv = serverIndexBegin[*itRank][1];
2028          const int ni_srv = serverDimensionSizes[*itRank][0];
2029          const int nj_srv = serverDimensionSizes[*itRank][1];
2030
2031          msgs.push_back(CMessage());
2032          CMessage& msg = msgs.back();
2033          msg << this->getId() ;
2034          msg << ni_srv << ibegin_srv << nj_srv << jbegin_srv;
2035          msg << global_zoom_ni.getValue() << global_zoom_ibegin.getValue() << global_zoom_nj.getValue() << global_zoom_jbegin.getValue();       
2036          msg << isCompressible_;
2037
2038          event.push(*itRank,1,msg);
2039        }
2040        contextClientTmp->sendEvent(event);
2041      }
2042      else contextClientTmp->sendEvent(event);
2043    }
2044  }
2045
2046  /*!
2047    Send mask index from client to connected(s)
2048  */
2049  void CDomain::sendMask()
2050  {
2051    int ns, n, i, j, ind, nv, idx;
2052    CContext* context = CContext::getCurrent();
2053
2054    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
2055    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
2056    for (int p = 0; p < nbSrvPools; ++p)
2057    {
2058      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
2059
2060      // send area for each connected server
2061      CEventClient eventMask(getType(), EVENT_ID_MASK);
2062
2063      list<CMessage> list_msgsMask;
2064      list<CArray<bool,1> > list_mask;
2065
2066      boost::unordered_map<int, vector<size_t> >::const_iterator it, iteMap;
2067      iteMap = indSrv_.end();
2068      for (int k = 0; k < connectedServerRank_.size(); ++k)
2069      {
2070        int nbData = 0;
2071        int rank = connectedServerRank_[k];
2072        it = indSrv_.find(rank);
2073        if (iteMap != it)
2074          nbData = it->second.size();
2075        list_mask.push_back(CArray<bool,1>(nbData));
2076
2077        const std::vector<size_t>& temp = it->second;
2078        for (n = 0; n < nbData; ++n)
2079        {
2080          idx = static_cast<int>(it->second[n]);
2081          list_mask.back()(n) = mask_1d(globalLocalIndexMap_[idx]);
2082        }
2083
2084        list_msgsMask.push_back(CMessage());
2085        list_msgsMask.back() << this->getId() << list_mask.back();
2086        eventMask.push(rank, nbConnectedClients_[rank], list_msgsMask.back());
2087      }
2088      client->sendEvent(eventMask);
2089    }
2090  }
2091
2092  /*!
2093    Send area from client to connected client(s)
2094  */
2095  void CDomain::sendArea()
2096  {
2097    if (!hasArea) return;
2098
2099    int ns, n, i, j, ind, nv, idx;
2100    CContext* context = CContext::getCurrent();
2101
2102    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
2103    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
2104    for (int p = 0; p < nbSrvPools; ++p)
2105    {
2106      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
2107
2108      // send area for each connected server
2109      CEventClient eventArea(getType(), EVENT_ID_AREA);
2110
2111      list<CMessage> list_msgsArea;
2112      list<CArray<double,1> > list_area;
2113
2114      boost::unordered_map<int, vector<size_t> >::const_iterator it, iteMap;
2115      iteMap = indSrv_.end();
2116      for (int k = 0; k < connectedServerRank_.size(); ++k)
2117      {
2118        int nbData = 0;
2119        int rank = connectedServerRank_[k];
2120        it = indSrv_.find(rank);
2121        if (iteMap != it)
2122          nbData = it->second.size();
2123        list_area.push_back(CArray<double,1>(nbData));
2124
2125        const std::vector<size_t>& temp = it->second;
2126        for (n = 0; n < nbData; ++n)
2127        {
2128          idx = static_cast<int>(it->second[n]);
2129          list_area.back()(n) = areavalue(globalLocalIndexMap_[idx]);
2130        }
2131
2132        list_msgsArea.push_back(CMessage());
2133        list_msgsArea.back() << this->getId() << hasArea;
2134        list_msgsArea.back() << list_area.back();
2135        eventArea.push(rank, nbConnectedClients_[rank], list_msgsArea.back());
2136      }
2137      client->sendEvent(eventArea);
2138    }
2139  }
2140
2141  /*!
2142    Send longitude and latitude from client to servers
2143    Each client send long and lat information to corresponding connected server(s).
2144    Because longitude and latitude are optional, this function only called if latitude and longitude exist
2145  */
2146  void CDomain::sendLonLat()
2147  {
2148    if (!hasLonLat) return;
2149
2150    int ns, n, i, j, ind, nv, idx;
2151    CContext* context = CContext::getCurrent();
2152
2153    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
2154    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
2155    for (int p = 0; p < nbSrvPools; ++p)
2156    {
2157      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
2158
2159      // send lon lat for each connected server
2160      CEventClient eventLon(getType(), EVENT_ID_LON);
2161      CEventClient eventLat(getType(), EVENT_ID_LAT);
2162
2163      list<CMessage> list_msgsLon, list_msgsLat;
2164      list<CArray<double,1> > list_lon, list_lat;
2165      list<CArray<double,2> > list_boundslon, list_boundslat;
2166
2167      boost::unordered_map<int, vector<size_t> >::const_iterator it, iteMap;
2168      iteMap = indSrv_.end();
2169      for (int k = 0; k < connectedServerRank_.size(); ++k)
2170      {
2171        int nbData = 0;
2172        int rank = connectedServerRank_[k];
2173        it = indSrv_.find(rank);
2174        if (iteMap != it)
2175          nbData = it->second.size();
2176
2177        list_lon.push_back(CArray<double,1>(nbData));
2178        list_lat.push_back(CArray<double,1>(nbData));
2179
2180        if (hasBounds)
2181        {
2182          list_boundslon.push_back(CArray<double,2>(nvertex, nbData));
2183          list_boundslat.push_back(CArray<double,2>(nvertex, nbData));
2184        }
2185
2186        CArray<double,1>& lon = list_lon.back();
2187        CArray<double,1>& lat = list_lat.back();
2188        const std::vector<size_t>& temp = it->second;
2189        for (n = 0; n < nbData; ++n)
2190        {
2191          idx = static_cast<int>(it->second[n]);
2192          int localInd = globalLocalIndexMap_[idx];
2193          lon(n) = lonvalue(localInd);
2194          lat(n) = latvalue(localInd);
2195
2196          if (hasBounds)
2197          {
2198            CArray<double,2>& boundslon = list_boundslon.back();
2199            CArray<double,2>& boundslat = list_boundslat.back();
2200
2201            for (nv = 0; nv < nvertex; ++nv)
2202            {
2203              boundslon(nv, n) = bounds_lonvalue(nv, localInd);
2204              boundslat(nv, n) = bounds_latvalue(nv, localInd);
2205            }
2206          }
2207        }
2208
2209        list_msgsLon.push_back(CMessage());
2210        list_msgsLat.push_back(CMessage());
2211
2212        list_msgsLon.back() << this->getId() << hasLonLat;
2213        if (hasLonLat) 
2214          list_msgsLon.back() << list_lon.back();
2215        list_msgsLon.back()  << hasBounds;
2216        if (hasBounds)
2217        {
2218          list_msgsLon.back() << list_boundslon.back();
2219        }
2220
2221        list_msgsLat.back() << this->getId() << hasLonLat;
2222        if (hasLonLat)
2223          list_msgsLat.back() << list_lat.back();
2224        list_msgsLat.back() << hasBounds;
2225        if (hasBounds)
2226        {         
2227          list_msgsLat.back() << list_boundslat.back();
2228        }
2229
2230        eventLon.push(rank, nbConnectedClients_[rank], list_msgsLon.back());
2231        eventLat.push(rank, nbConnectedClients_[rank], list_msgsLat.back());
2232      }
2233      client->sendEvent(eventLon);
2234      client->sendEvent(eventLat);
2235    }
2236  }
2237
2238  /*!
2239    Send data index to corresponding connected clients.
2240    Data index can be compressed however, we always send decompressed data index
2241    and they will be compressed on receiving.
2242  */
2243  void CDomain::sendDataIndex()
2244  {
2245    int ns, n, i, j, ind, nv, idx;
2246    CContext* context = CContext::getCurrent();
2247
2248    // int nbSrvPools = (context->hasServer) ? context->clientPrimServer.size() : 1;
2249    int nbSrvPools = (context->hasServer) ? (context->hasClient ? context->clientPrimServer.size() : 0) : 1;
2250    for (int p = 0; p < nbSrvPools; ++p)
2251    {
2252      CContextClient* client = (0 != context->clientPrimServer.size()) ? context->clientPrimServer[p] : context->client;
2253
2254      // send area for each connected server
2255      CEventClient eventDataIndex(getType(), EVENT_ID_DATA_INDEX);
2256
2257      list<CMessage> list_msgsDataIndex;
2258      list<CArray<int,1> > list_data_i_index, list_data_j_index;
2259
2260      int nbIndex = i_index.numElements();
2261      int niByIndex = max(i_index) - min(i_index) + 1;
2262      int njByIndex = max(j_index) - min(j_index) + 1; 
2263      int dataIindexBound = (1 == data_dim) ? (niByIndex * njByIndex) : niByIndex;
2264      int dataJindexBound = (1 == data_dim) ? (niByIndex * njByIndex) : njByIndex;
2265
2266     
2267      CArray<int,1> dataIIndex(nbIndex), dataJIndex(nbIndex);
2268      dataIIndex = -1; 
2269      dataJIndex = -1;
2270      ind = 0;
2271
2272      for (idx = 0; idx < data_i_index.numElements(); ++idx)
2273      {
2274        int dataIidx = data_i_index(idx) + data_ibegin;
2275        int dataJidx = data_j_index(idx) + data_jbegin;
2276        if ((0 <= dataIidx) && (dataIidx < dataIindexBound) &&
2277            (0 <= dataJidx) && (dataJidx < dataJindexBound))
2278        {
2279          dataIIndex((1 == data_dim) ? dataIidx : dataJidx * ni + dataIidx) = 1; //i_index(dataIidx);//dataIidx;
2280          dataJIndex((1 == data_dim) ? dataIidx : dataJidx * ni + dataIidx) = 1; //j_index(dataJidx);//         
2281        }
2282      }
2283
2284      boost::unordered_map<int, vector<size_t> >::const_iterator it, iteMap;
2285      iteMap = indSrv_.end();
2286      for (int k = 0; k < connectedServerRank_.size(); ++k)
2287      {
2288        int nbData = 0;
2289        int rank = connectedServerRank_[k];
2290        it = indSrv_.find(rank);
2291        if (iteMap != it)
2292          nbData = it->second.size();
2293        list_data_i_index.push_back(CArray<int,1>(nbData));
2294        list_data_j_index.push_back(CArray<int,1>(nbData));
2295
2296        const std::vector<size_t>& temp = it->second;
2297        for (n = 0; n < nbData; ++n)
2298        {
2299          idx = static_cast<int>(it->second[n]);
2300          i = globalLocalIndexMap_[idx];
2301          list_data_i_index.back()(n) = dataIIndex(i);
2302          list_data_j_index.back()(n) = dataJIndex(i);
2303        }
2304
2305        list_msgsDataIndex.push_back(CMessage());
2306        list_msgsDataIndex.back() << this->getId();
2307        list_msgsDataIndex.back() << list_data_i_index.back() << list_data_j_index.back();
2308        eventDataIndex.push(rank, nbConnectedClients_[rank], list_msgsDataIndex.back());
2309      }
2310      client->sendEvent(eventDataIndex);
2311    }
2312  }
2313 
2314  bool CDomain::dispatchEvent(CEventServer& event)
2315  {
2316    if (SuperClass::dispatchEvent(event)) return true;
2317    else
2318    {
2319      switch(event.type)
2320      {
2321        case EVENT_ID_SERVER_ATTRIBUT:
2322          recvDistributionAttributes(event);
2323          return true;
2324          break;
2325        case EVENT_ID_INDEX:
2326          recvIndex(event);
2327          return true;
2328          break;
2329        case EVENT_ID_INDEX_ZOOM:
2330          recvIndexZoom(event);
2331          return true;
2332          break;
2333        case EVENT_ID_MASK:
2334          recvMask(event);
2335          return true;
2336          break;
2337        case EVENT_ID_LON:
2338          recvLon(event);
2339          return true;
2340          break;
2341        case EVENT_ID_LAT:
2342          recvLat(event);
2343          return true;
2344          break;
2345        case EVENT_ID_AREA:
2346          recvArea(event);
2347          return true;
2348          break; 
2349        case EVENT_ID_DATA_INDEX:
2350          recvDataIndex(event);
2351          return true;
2352          break;
2353        default:
2354          ERROR("bool CDomain::dispatchEvent(CEventServer& event)",
2355                << "Unknown Event");
2356          return false;
2357       }
2358    }
2359  }
2360
2361  /*!
2362    Receive index event from clients(s)
2363    \param[in] event event contain info about rank and associated index
2364  */
2365  void CDomain::recvIndex(CEventServer& event)
2366  {
2367    string domainId;
2368    std::map<int, CBufferIn*> rankBuffers;
2369
2370    list<CEventServer::SSubEvent>::iterator it;
2371    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2372    {     
2373      CBufferIn* buffer = it->buffer;
2374      *buffer >> domainId;
2375      rankBuffers[it->rank] = buffer;       
2376    }
2377    get(domainId)->recvIndex(rankBuffers);
2378
2379    // if (domain->isCompressible_)
2380    // {
2381    //   std::sort(domain->indexesToWrite.begin(), domain->indexesToWrite.end());
2382
2383    //   CContextServer* server = CContext::getCurrent()->server;
2384    //   domain->numberWrittenIndexes_ = domain->indexesToWrite.size();
2385    //   MPI_Allreduce(&domain->numberWrittenIndexes_, &domain->totalNumberWrittenIndexes_, 1, MPI_INT, MPI_SUM, server->intraComm);
2386    //   MPI_Scan(&domain->numberWrittenIndexes_, &domain->offsetWrittenIndexes_, 1, MPI_INT, MPI_SUM, server->intraComm);
2387    //   domain->offsetWrittenIndexes_ -= domain->numberWrittenIndexes_;
2388    // }
2389  }
2390
2391  /*!
2392    Receive index information from client(s)
2393    \param[in] rankBuffers rank of sending client and the corresponding receive buffer 
2394  */
2395  void CDomain::recvIndex(std::map<int, CBufferIn*>& rankBuffers)
2396  {
2397    int nbReceived = rankBuffers.size(), i, ind, index, type_int;
2398    recvClientRanks_.resize(nbReceived);       
2399
2400    std::map<int, CBufferIn*>::iterator it = rankBuffers.begin(), ite = rankBuffers.end();
2401    ind = 0;
2402    for (ind = 0; it != ite; ++it, ++ind)
2403    {       
2404       recvClientRanks_[ind] = it->first;
2405       CBufferIn& buffer = *(it->second);
2406       buffer >> type_int >> isCurvilinear >> indGlob_[it->first]; 
2407       type.setValue((type_attr::t_enum)type_int); // probleme des type enum avec les buffers : ToFix
2408    }
2409    int nbIndGlob = 0;
2410    for (i = 0; i < nbReceived; ++i)
2411    {
2412      nbIndGlob += indGlob_[recvClientRanks_[i]].numElements();
2413    }
2414   
2415    i_index.resize(nbIndGlob);
2416    j_index.resize(nbIndGlob);
2417
2418    globalLocalIndexMap_.rehash(std::ceil(nbIndGlob/globalLocalIndexMap_.max_load_factor()));
2419    nbIndGlob = 0;
2420    for (i = 0; i < nbReceived; ++i)
2421    {
2422      CArray<int,1>& tmp = indGlob_[recvClientRanks_[i]];
2423      for (ind = 0; ind < tmp.numElements(); ++ind)
2424      {
2425         index = tmp(ind);
2426         i_index(nbIndGlob) = index % ni_glo;
2427         j_index(nbIndGlob) = index / ni_glo;
2428         globalLocalIndexMap_[index] = nbIndGlob; 
2429         ++nbIndGlob;
2430      } 
2431    }   
2432   
2433    // {
2434    //   CContextServer* server = CContext::getCurrent()->server;
2435    //   count_write_index_.resize(2);
2436    //   start_write_index_.resize(2);
2437    //   local_write_size_.resize(2);
2438    //   global_write_size_.resize(2);
2439    //   if ((this->type) == CDomain::type_attr::unstructured)
2440    //   {
2441    //     count_write_index_[0] = zoom_i_index.numElements();
2442    //     count_write_index_[1] = 0;
2443    //   }
2444    //   else
2445    //   {
2446    //     int ni_zoom = zoom_i_index.numElements(), idx, nbIZoom = 0, nbJZoom = 0;
2447    //     for (idx =0; idx < ni_zoom; ++idx)
2448    //     {
2449    //        if ((ibegin <= zoom_i_index(idx)) && (zoom_i_index(idx) < ibegin+ni) && (nbIZoom < ni))
2450    //         ++nbIZoom;
2451    //        if ((jbegin <= zoom_j_index(idx)) && (zoom_j_index(idx) < jbegin+nj) && (nbJZoom < nj))
2452    //         ++nbJZoom;
2453    //     }
2454    //     count_write_index_[0] = nbIZoom;
2455    //     count_write_index_[1] = nbJZoom;
2456
2457    //     // Reoder the zoom_index
2458    //     for (int j = 0; j < nbJZoom; ++j)
2459    //       for (int i = 0; i < nbIZoom; ++i)
2460    //       {
2461    //         idx = nbIZoom * j + i;
2462    //         if (idx < ni_zoom)
2463    //         {
2464    //           zoom_i_index(idx) = ibegin + i;
2465    //           zoom_j_index(idx) = jbegin + j;
2466    //         }
2467    //       } 
2468
2469    //     // Reorder the global index
2470    //     for (int j = 0; j < nj; ++j)
2471    //       for (int i = 0; i < ni; ++i)
2472    //       {
2473    //         idx = ni * j + i;
2474    //         if (idx < nbIndGlob)
2475    //         {
2476    //           i_index(idx) = ibegin + i;
2477    //           j_index(idx) = jbegin + j;
2478    //         }
2479    //       }         
2480    //   }
2481
2482           
2483    //   MPI_Scan(&count_write_index_[0], &start_write_index_[0], 2, MPI_INT, MPI_SUM, server->intraComm);     
2484    //   start_write_index_[0] = 0;
2485    //   start_write_index_[1] -= count_write_index_[1];
2486    //   local_write_size_[0] = count_write_index_[0];
2487    //   local_write_size_[1] = count_write_index_[1];
2488    //   MPI_Allreduce(&count_write_index_[0], &global_write_size_[0], 2, MPI_INT, MPI_SUM, server->intraComm);
2489    //   global_write_size_[0] = count_write_index_[0];
2490    //   global_write_size_[1] = (global_write_size_[1] > nj_glo) ? nj_glo : global_write_size_[1];
2491         
2492
2493    // }
2494
2495    // int type_int;
2496    // buffer >> type_int >> isCurvilinear >> indiSrv[rank] >> indjSrv[rank];
2497    // type.setValue((type_attr::t_enum)type_int); // probleme des type enum avec les buffers : ToFix
2498
2499    // if (isCompressible_)
2500    // {
2501    //   CArray<int, 1> writtenIndexes;
2502    //   buffer >> writtenIndexes;
2503    //   indexesToWrite.reserve(indexesToWrite.size() + writtenIndexes.numElements());
2504    //   for (int i = 0; i < writtenIndexes.numElements(); ++i)
2505    //     indexesToWrite.push_back(writtenIndexes(i));
2506    // }
2507  }
2508
2509  /*!
2510    Receive index event from clients(s)
2511    \param[in] event event contain info about rank and associated index
2512  */
2513  void CDomain::recvIndexZoom(CEventServer& event)
2514  {
2515    string domainId;
2516    std::map<int, CBufferIn*> rankBuffers;
2517
2518    list<CEventServer::SSubEvent>::iterator it;
2519    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2520    {     
2521      CBufferIn* buffer = it->buffer;
2522      *buffer >> domainId;
2523      rankBuffers[it->rank] = buffer;       
2524    }
2525    get(domainId)->recvIndexZoom(rankBuffers);
2526  }
2527
2528  /*!
2529    Receive index information from client(s)
2530    \param[in] rankBuffers rank of sending client and the corresponding receive buffer 
2531  */
2532  void CDomain::recvIndexZoom(std::map<int, CBufferIn*>& rankBuffers)
2533  {
2534    int nbReceived = rankBuffers.size(), i, ind, index, type_int;
2535    recvClientZoomRanks_.resize(nbReceived);   
2536    int ni_zoom_tmp, ibegin_zoom_tmp, nj_zoom_tmp, jbegin_zoom_tmp;
2537
2538    std::map<int, CBufferIn*>::iterator it = rankBuffers.begin(), ite = rankBuffers.end();
2539    ind = 0;
2540    for (ind = 0; it != ite; ++it, ++ind)
2541    {       
2542       recvClientZoomRanks_[ind] = it->first;
2543       CBufferIn& buffer = *(it->second);
2544       buffer >> indGlobZoom_[it->first] >> doZoomByIndex_;       
2545    }
2546
2547    int nbZoomInd = 0;
2548    for (i = 0; i < nbReceived; ++i)
2549    {
2550      nbZoomInd += indGlobZoom_[recvClientZoomRanks_[i]].numElements();
2551    }
2552
2553    if (doZoomByIndex_)
2554    {
2555      zoom_i_index.resize(nbZoomInd);
2556      zoom_j_index.resize(nbZoomInd);
2557     
2558      nbZoomInd = 0;
2559      for (i = 0; i < nbReceived; ++i)
2560      {
2561        CArray<int,1>& tmp = indGlobZoom_[recvClientRanks_[i]];
2562        for (ind = 0; ind < tmp.numElements(); ++ind)
2563        {
2564           index = tmp(ind);
2565           zoom_i_index(nbZoomInd) = index % ni_glo;
2566           zoom_j_index(nbZoomInd) = index / ni_glo;
2567           ++nbZoomInd;
2568        } 
2569      }     
2570    }
2571    else 
2572    {
2573    }
2574
2575    // globalLocalIndexZoomMap_.rehash(std::ceil(nbZoomInd/globalLocalIndexZoomMap_.max_load_factor()));
2576    // nbZoomInd = 0;
2577    // for (int j = 0; j < zoom_nj; ++j)
2578    //   for (int i = 0; i < zoom_ni; ++i)
2579    //   {
2580    //     globalLocalIndexZoomMap_[(i + zoom_ibegin) + (j + zoom_jbegin) * ni_glo] = nbZoomInd;
2581    //     ++nbZoomInd;
2582    //   }
2583  }
2584
2585  /*!
2586    Receive attributes event from clients(s)
2587    \param[in] event event contain info about rank and associated attributes
2588  */
2589  void CDomain::recvDistributionAttributes(CEventServer& event)
2590  {
2591    CBufferIn* buffer=event.subEvents.begin()->buffer;
2592    string domainId ;
2593    *buffer>>domainId ;
2594    get(domainId)->recvDistributionAttributes(*buffer);
2595  }
2596
2597  /*!
2598    Receive attributes from client(s): zoom info and begin and n of each server
2599    \param[in] rank rank of client source
2600    \param[in] buffer message containing attributes info
2601  */
2602  void CDomain::recvDistributionAttributes(CBufferIn& buffer)
2603  {
2604    int ni_tmp, ibegin_tmp, nj_tmp, jbegin_tmp;
2605    int global_zoom_ni_tmp, global_zoom_ibegin_tmp, global_zoom_nj_tmp, global_zoom_jbegin_tmp;
2606    buffer >> ni_tmp >> ibegin_tmp >> nj_tmp >> jbegin_tmp
2607           >> global_zoom_ni_tmp >> global_zoom_ibegin_tmp >> global_zoom_nj_tmp >> global_zoom_jbegin_tmp           
2608           >> isCompressible_;
2609    ni.setValue(ni_tmp);
2610    ibegin.setValue(ibegin_tmp);
2611    nj.setValue(nj_tmp);
2612    jbegin.setValue(jbegin_tmp);
2613
2614    global_zoom_ni.setValue(global_zoom_ni_tmp);
2615    global_zoom_ibegin.setValue(global_zoom_ibegin_tmp);
2616    global_zoom_nj.setValue(global_zoom_nj_tmp);
2617    global_zoom_jbegin.setValue(global_zoom_jbegin_tmp);
2618
2619    int iend = ibegin + ni  - 1;
2620    int jend = jbegin + nj  - 1;
2621    int zoom_iend_glob = global_zoom_ibegin + global_zoom_ni - 1;
2622    int zoom_jend_glob = global_zoom_jbegin + global_zoom_nj - 1;
2623
2624    zoom_ibegin.setValue(global_zoom_ibegin > ibegin ? global_zoom_ibegin : ibegin);
2625    int zoom_iend = zoom_iend_glob < iend ? zoom_iend_glob : iend ;
2626    zoom_ni.setValue(zoom_iend-zoom_ibegin+1);
2627
2628    zoom_jbegin.setValue(global_zoom_jbegin > jbegin ? global_zoom_jbegin : jbegin);
2629    int zoom_jend = zoom_jend_glob < jend ? zoom_jend_glob : jend ;
2630    zoom_nj.setValue(zoom_jend-zoom_jbegin+1);
2631
2632    if (zoom_ni<=0 || zoom_nj<=0)
2633    {
2634      zoom_ibegin=0 ; zoom_iend=0 ; zoom_ni=0 ;
2635      zoom_jbegin=0 ; zoom_jend=0 ; zoom_nj=0 ;
2636    }
2637
2638  }
2639
2640  /*!
2641    Receive area event from clients(s)
2642    \param[in] event event contain info about rank and associated area
2643  */
2644  void CDomain::recvMask(CEventServer& event)
2645  {
2646    string domainId;
2647    std::map<int, CBufferIn*> rankBuffers;
2648
2649    list<CEventServer::SSubEvent>::iterator it;
2650    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2651    {     
2652      CBufferIn* buffer = it->buffer;
2653      *buffer >> domainId;
2654      rankBuffers[it->rank] = buffer;     
2655    }
2656    get(domainId)->recvMask(rankBuffers);
2657  }
2658
2659
2660  /*!
2661    Receive mask information from client(s)
2662    \param[in] rankBuffers rank of sending client and the corresponding receive buffer 
2663  */
2664  void CDomain::recvMask(std::map<int, CBufferIn*>& rankBuffers)
2665  {
2666    int nbReceived = rankBuffers.size(), i, ind, index;
2667    if (nbReceived != recvClientRanks_.size())
2668      ERROR("void CDomain::recvArea(std::map<int, CBufferIn*>& rankBuffers)",
2669           << "The number of sending clients is not correct."
2670           << "Expected number: " << recvClientRanks_.size() << " but received " << nbReceived);
2671
2672    vector<CArray<bool,1> > recvMaskValue(nbReceived);     
2673    for (i = 0; i < recvClientRanks_.size(); ++i)
2674    {
2675      int rank = recvClientRanks_[i];
2676      CBufferIn& buffer = *(rankBuffers[rank]);     
2677      buffer >> recvMaskValue[i];
2678    }
2679
2680    int nbMaskInd = 0;
2681    for (i = 0; i < nbReceived; ++i)
2682    {
2683      nbMaskInd += recvMaskValue[i].numElements();
2684    }
2685 
2686    mask_1d.resize(nbMaskInd);
2687    nbMaskInd = 0;
2688    for (i = 0; i < nbReceived; ++i)
2689    {
2690      CArray<bool,1>& tmp = recvMaskValue[i];
2691      for (ind = 0; ind < tmp.numElements(); ++ind)
2692      {
2693        mask_1d(nbMaskInd) = tmp(ind);     
2694        ++nbMaskInd;
2695      }
2696    }   
2697  }
2698
2699  /*!
2700    Receive longitude event from clients(s)
2701    \param[in] event event contain info about rank and associated longitude
2702  */
2703  void CDomain::recvLon(CEventServer& event)
2704  {
2705    string domainId;
2706    std::map<int, CBufferIn*> rankBuffers;
2707
2708    list<CEventServer::SSubEvent>::iterator it;
2709    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2710    {     
2711      CBufferIn* buffer = it->buffer;
2712      *buffer >> domainId;
2713      rankBuffers[it->rank] = buffer;       
2714    }
2715    get(domainId)->recvLon(rankBuffers);
2716  }
2717
2718  /*!
2719    Receive longitude information from client(s)
2720    \param[in] rankBuffers rank of sending client and the corresponding receive buffer 
2721  */
2722  void CDomain::recvLon(std::map<int, CBufferIn*>& rankBuffers)
2723  {
2724    int nbReceived = rankBuffers.size(), i, ind, index, iindex, jindex, lInd;
2725    if (nbReceived != recvClientRanks_.size())
2726      ERROR("void CDomain::recvLon(std::map<int, CBufferIn*>& rankBuffers)",
2727           << "The number of sending clients is not correct."
2728           << "Expected number: " << recvClientRanks_.size() << " but received " << nbReceived);
2729
2730    vector<CArray<double,1> > recvLonValue(nbReceived);
2731    vector<CArray<double,2> > recvBoundsLonValue(nbReceived);   
2732    for (i = 0; i < recvClientRanks_.size(); ++i)
2733    {
2734      int rank = recvClientRanks_[i];
2735      CBufferIn& buffer = *(rankBuffers[rank]);
2736      buffer >> hasLonLat;
2737      if (hasLonLat)
2738        buffer >> recvLonValue[i];
2739      buffer >> hasBounds;
2740      if (hasBounds)
2741        buffer >> recvBoundsLonValue[i];
2742    }
2743
2744    if (hasLonLat)
2745    {
2746      int nbLonInd = 0;
2747      for (i = 0; i < nbReceived; ++i)
2748      {
2749        nbLonInd += recvLonValue[i].numElements();
2750      }
2751   
2752      if (nbLonInd != globalLocalIndexMap_.size())
2753        info (0) << "Something wrong with longitude index "<< std::endl;
2754
2755      lonvalue.resize(nbLonInd);
2756      if (hasBounds)
2757      {
2758        bounds_lonvalue.resize(nvertex,nbLonInd);
2759        bounds_lonvalue = 0.;
2760      }
2761
2762      nbLonInd = 0;
2763      for (i = 0; i < nbReceived; ++i)
2764      {
2765        CArray<int,1>& tmpInd = indGlob_[recvClientRanks_[i]];
2766        CArray<double,1>& tmp = recvLonValue[i];
2767        for (ind = 0; ind < tmp.numElements(); ++ind)
2768        {
2769          lInd = globalLocalIndexMap_[size_t(tmpInd(ind))];
2770          lonvalue(lInd) = tmp(ind); 
2771           if (hasBounds)
2772           {         
2773            for (int nv = 0; nv < nvertex; ++nv)
2774              bounds_lonvalue(nv, lInd) = recvBoundsLonValue[i](nv, ind);
2775           }                 
2776        }
2777      }       
2778    }
2779  }
2780
2781  /*!
2782    Receive latitude event from clients(s)
2783    \param[in] event event contain info about rank and associated latitude
2784  */
2785  void CDomain::recvLat(CEventServer& event)
2786  {
2787    string domainId;
2788    std::map<int, CBufferIn*> rankBuffers;
2789
2790    list<CEventServer::SSubEvent>::iterator it;
2791    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2792    {     
2793      CBufferIn* buffer = it->buffer;
2794      *buffer >> domainId;
2795      rankBuffers[it->rank] = buffer;   
2796    }
2797    get(domainId)->recvLat(rankBuffers);
2798  }
2799
2800  /*!
2801    Receive latitude information from client(s)
2802    \param[in] rankBuffers rank of sending client and the corresponding receive buffer 
2803  */
2804  void CDomain::recvLat(std::map<int, CBufferIn*>& rankBuffers)
2805  {
2806    int nbReceived = rankBuffers.size(), i, ind, index, iindex, jindex, lInd;
2807    if (nbReceived != recvClientRanks_.size())
2808      ERROR("void CDomain::recvLat(std::map<int, CBufferIn*>& rankBuffers)",
2809           << "The number of sending clients is not correct."
2810           << "Expected number: " << recvClientRanks_.size() << " but received " << nbReceived);
2811
2812    vector<CArray<double,1> > recvLatValue(nbReceived);
2813    vector<CArray<double,2> > recvBoundsLatValue(nbReceived);   
2814    for (i = 0; i < recvClientRanks_.size(); ++i)
2815    {
2816      int rank = recvClientRanks_[i];
2817      CBufferIn& buffer = *(rankBuffers[rank]);
2818      buffer >> hasLonLat;
2819      if (hasLonLat)
2820        buffer >> recvLatValue[i];
2821      buffer >> hasBounds;
2822      if (hasBounds)
2823        buffer >> recvBoundsLatValue[i];
2824    }
2825
2826    if (hasLonLat)
2827    {
2828      int nbLatInd = 0;
2829      for (i = 0; i < nbReceived; ++i)
2830      {
2831        nbLatInd += recvLatValue[i].numElements();
2832      }
2833   
2834      if (nbLatInd != globalLocalIndexMap_.size())
2835        info (0) << "Something wrong with latitude index "<< std::endl;
2836
2837      latvalue.resize(nbLatInd);
2838      if (hasBounds)
2839      {
2840        bounds_latvalue.resize(nvertex,nbLatInd);
2841        bounds_latvalue = 0. ;
2842      }
2843
2844      nbLatInd = 0;
2845      for (i = 0; i < nbReceived; ++i)
2846      {
2847        CArray<int,1>& tmpInd = indGlob_[recvClientRanks_[i]];
2848        CArray<double,1>& tmp = recvLatValue[i];
2849        for (ind = 0; ind < tmp.numElements(); ++ind)
2850        {
2851          lInd = globalLocalIndexMap_[size_t(tmpInd(ind))];
2852          latvalue(lInd) = tmp(ind);   
2853           if (hasBounds)
2854           {
2855            CArray<double,2>& boundslat = recvBoundsLatValue[i];
2856            for (int nv = 0; nv < nvertex; ++nv)
2857              bounds_latvalue(nv, lInd) = boundslat(nv, ind);
2858           }   
2859          ++nbLatInd;
2860        }
2861      }       
2862    }
2863  }
2864
2865  /*!
2866    Receive area event from clients(s)
2867    \param[in] event event contain info about rank and associated area
2868  */
2869  void CDomain::recvArea(CEventServer& event)
2870  {
2871    string domainId;
2872    std::map<int, CBufferIn*> rankBuffers;
2873
2874    list<CEventServer::SSubEvent>::iterator it;
2875    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2876    {     
2877      CBufferIn* buffer = it->buffer;
2878      *buffer >> domainId;
2879      rankBuffers[it->rank] = buffer;     
2880    }
2881    get(domainId)->recvArea(rankBuffers);
2882  }
2883
2884
2885  /*!
2886    Receive area information from client(s)
2887    \param[in] rankBuffers rank of sending client and the corresponding receive buffer     
2888  */
2889  void CDomain::recvArea(std::map<int, CBufferIn*>& rankBuffers)
2890  {
2891    int nbReceived = rankBuffers.size(), i, ind, index, lInd;
2892    if (nbReceived != recvClientRanks_.size())
2893      ERROR("void CDomain::recvArea(std::map<int, CBufferIn*>& rankBuffers)",
2894           << "The number of sending clients is not correct."
2895           << "Expected number: " << recvClientRanks_.size() << " but received " << nbReceived);
2896
2897    vector<CArray<double,1> > recvAreaValue(nbReceived);     
2898    for (i = 0; i < recvClientRanks_.size(); ++i)
2899    {
2900      int rank = recvClientRanks_[i];
2901      CBufferIn& buffer = *(rankBuffers[rank]);     
2902      buffer >> hasArea;
2903      if (hasArea)
2904        buffer >> recvAreaValue[i];
2905    }
2906
2907    int nbAreaInd = 0;
2908    for (i = 0; i < nbReceived; ++i)
2909    {     
2910      nbAreaInd += recvAreaValue[i].numElements();
2911    }
2912 
2913    if (nbAreaInd != globalLocalIndexMap_.size())
2914      info (0) << "Something wrong with latitude index "<< std::endl;
2915
2916    if (hasArea)
2917    {
2918      areavalue.resize(nbAreaInd);
2919      nbAreaInd = 0;     
2920      for (i = 0; i < nbReceived; ++i)
2921      {
2922        CArray<int,1>& tmpInd = indGlob_[recvClientRanks_[i]];
2923        CArray<double,1>& tmp = recvAreaValue[i];
2924        for (ind = 0; ind < tmp.numElements(); ++ind)
2925        {
2926          lInd = globalLocalIndexMap_[size_t(tmpInd(ind))];
2927          areavalue(lInd) = tmp(ind);         
2928        }
2929      }
2930     
2931    }
2932  }
2933
2934  /*!
2935    Receive data index event from clients(s)
2936    \param[in] event event contain info about rank and associated index
2937  */
2938  void CDomain::recvDataIndex(CEventServer& event)
2939  {
2940    string domainId;
2941    std::map<int, CBufferIn*> rankBuffers;
2942
2943    list<CEventServer::SSubEvent>::iterator it;
2944    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
2945    {     
2946      CBufferIn* buffer = it->buffer;
2947      *buffer >> domainId;
2948      rankBuffers[it->rank] = buffer;       
2949    }
2950    get(domainId)->recvDataIndex(rankBuffers);
2951
2952    // if (domain->isCompressible_)
2953    // {
2954    //   std::sort(domain->indexesToWrite.begin(), domain->indexesToWrite.end());
2955
2956    //   CContextServer* server = CContext::getCurrent()->server;
2957    //   domain->numberWrittenIndexes_ = domain->indexesToWrite.size();
2958    //   MPI_Allreduce(&domain->numberWrittenIndexes_, &domain->totalNumberWrittenIndexes_, 1, MPI_INT, MPI_SUM, server->intraComm);
2959    //   MPI_Scan(&domain->numberWrittenIndexes_, &domain->offsetWrittenIndexes_, 1, MPI_INT, MPI_SUM, server->intraComm);
2960    //   domain->offsetWrittenIndexes_ -= domain->numberWrittenIndexes_;
2961    // }
2962  }
2963
2964  /*!
2965    Receive data index information from client(s)
2966    A client receives data index from different clients to rebuild its own data index.
2967    Because the data index is local, to rebuild data index of received client, we should use global index along with.
2968    \param[in] rankBuffers rank of sending client and the corresponding receive buffer     
2969  */
2970  void CDomain::recvDataIndex(std::map<int, CBufferIn*>& rankBuffers)
2971  {
2972    int nbReceived = rankBuffers.size(), i, ind, index, indexI, indexJ, type_int, lInd;   
2973    if (nbReceived != recvClientRanks_.size())
2974      ERROR("void CDomain::recvDataIndex(std::map<int, CBufferIn*>& rankBuffers)",
2975           << "The number of sending clients is not correct."
2976           << "Expected number: " << recvClientRanks_.size() << " but received " << nbReceived);
2977
2978    vector<CArray<int,1> > recvDataIIndex(nbReceived),recvDataJIndex(nbReceived);     
2979    for (i = 0; i < recvClientRanks_.size(); ++i)
2980    {
2981      int rank = recvClientRanks_[i];
2982      CBufferIn& buffer = *(rankBuffers[rank]);
2983      buffer >> recvDataIIndex[i];
2984      buffer >> recvDataJIndex[i];
2985    }
2986   
2987    int nbIndex = i_index.numElements();
2988    CArray<int,1> dataIIndex(nbIndex), dataJIndex(nbIndex);
2989    dataIIndex = -1; dataJIndex = -1;
2990
2991    nbIndex = 0;
2992    for (i = 0; i < nbReceived; ++i)
2993    {     
2994      CArray<int,1>& tmpInd = indGlob_[recvClientRanks_[i]];
2995      CArray<int,1>& tmpI = recvDataIIndex[i];   
2996      CArray<int,1>& tmpJ = recvDataJIndex[i];     
2997      if ((tmpI.numElements() != tmpInd.numElements()) || (tmpJ.numElements() != tmpInd.numElements()))
2998          ERROR("void CDomain::recvDataIndex(std::map<int, CBufferIn*>& rankBuffers)",
2999             << "The number of global received index is not coherent with the number of received data index."
3000             << "Expected number of global index: " << tmpI.numElements() << " but received " << tmpInd.numElements());
3001
3002      for (ind = 0; ind < tmpI.numElements(); ++ind)
3003      {
3004         lInd = globalLocalIndexMap_[size_t(tmpInd(ind))];
3005         dataIIndex(lInd) = tmpI(ind);
3006         dataJIndex(lInd) = tmpJ(ind);         
3007      } 
3008    }
3009
3010    int nbCompressedData = 0; 
3011    for (ind = 0; ind < dataIIndex.numElements(); ++ind)
3012    {
3013       indexI = dataIIndex(ind); indexJ = dataJIndex(ind);
3014       if ((0 <= indexI) && (0 <= indexJ))
3015         ++nbCompressedData;
3016    }       
3017 
3018    data_i_index.resize(nbCompressedData);
3019    data_j_index.resize(nbCompressedData);
3020
3021    nbCompressedData = 0; 
3022    for (ind = 0; ind < dataIIndex.numElements(); ++ind)
3023    {
3024       indexI = dataIIndex(ind); indexJ = dataJIndex(ind);
3025       if ((0 <= indexI) && (0 <= indexJ))
3026       {
3027          data_i_index(nbCompressedData) = (1 == data_dim) ? ind : ind % ni;
3028          data_j_index(nbCompressedData) = (1 == data_dim) ? 0   : ind / ni;
3029         ++nbCompressedData;
3030       }
3031    }
3032
3033    // Reset data_ibegin, data_jbegin
3034    data_ibegin.setValue(0);
3035    data_jbegin.setValue(0);
3036  }
3037
3038  CTransformation<CDomain>* CDomain::addTransformation(ETranformationType transType, const StdString& id)
3039  {
3040    transformationMap_.push_back(std::make_pair(transType, CTransformation<CDomain>::createTransformation(transType,id)));
3041    return transformationMap_.back().second;
3042  }
3043
3044  /*!
3045    Check whether a domain has transformation
3046    \return true if domain has transformation
3047  */
3048  bool CDomain::hasTransformation()
3049  {
3050    return (!transformationMap_.empty());
3051  }
3052
3053  /*!
3054    Set transformation for current domain. It's the method to move transformation in hierarchy
3055    \param [in] domTrans transformation on domain
3056  */
3057  void CDomain::setTransformations(const TransMapTypes& domTrans)
3058  {
3059    transformationMap_ = domTrans;
3060  }
3061
3062  /*!
3063    Get all transformation current domain has
3064    \return all transformation
3065  */
3066  CDomain::TransMapTypes CDomain::getAllTransformations(void)
3067  {
3068    return transformationMap_;
3069  }
3070
3071  void CDomain::duplicateTransformation(CDomain* src)
3072  {
3073    if (src->hasTransformation())
3074    {
3075      this->setTransformations(src->getAllTransformations());
3076    }
3077  }
3078
3079  /*!
3080   * Go through the hierarchy to find the domain from which the transformations must be inherited
3081   */
3082  void CDomain::solveInheritanceTransformation()
3083  {
3084    if (hasTransformation() || !hasDirectDomainReference())
3085      return;
3086
3087    CDomain* domain = this;
3088    std::vector<CDomain*> refDomains;
3089    while (!domain->hasTransformation() && domain->hasDirectDomainReference())
3090    {
3091      refDomains.push_back(domain);
3092      domain = domain->getDirectDomainReference();
3093    }
3094
3095    if (domain->hasTransformation())
3096      for (size_t i = 0; i < refDomains.size(); ++i)
3097        refDomains[i]->setTransformations(domain->getAllTransformations());
3098  }
3099
3100  /*!
3101    Parse children nodes of a domain in xml file.
3102    Whenver there is a new transformation, its type and name should be added into this function
3103    \param node child node to process
3104  */
3105  void CDomain::parse(xml::CXMLNode & node)
3106  {
3107    SuperClass::parse(node);
3108
3109    if (node.goToChildElement())
3110    {
3111      StdString nodeElementName;
3112      do
3113      {
3114        StdString nodeId("");
3115        if (node.getAttributes().end() != node.getAttributes().find("id"))
3116        { nodeId = node.getAttributes()["id"]; }
3117
3118        nodeElementName = node.getElementName();
3119        std::map<StdString, ETranformationType>::const_iterator ite = transformationMapList_.end(), it;
3120        it = transformationMapList_.find(nodeElementName);
3121        if (ite != it)
3122        {
3123          transformationMap_.push_back(std::make_pair(it->second, CTransformation<CDomain>::createTransformation(it->second,
3124                                                                                                                nodeId,
3125                                                                                                                &node)));
3126        }
3127        else
3128        {
3129          ERROR("void CDomain::parse(xml::CXMLNode & node)",
3130                << "The transformation " << nodeElementName << " has not been supported yet.");
3131        }
3132      } while (node.goToNextElement()) ;
3133      node.goToParentElement();
3134    }
3135  }
3136   //----------------------------------------------------------------
3137
3138   DEFINE_REF_FUNC(Domain,domain)
3139
3140   ///---------------------------------------------------------------
3141
3142} // namespace xios
Note: See TracBrowser for help on using the repository browser.