source: XIOS/dev/XIOS_DEV_CMIP6/src/node/domain.cpp @ 1431

Last change on this file since 1431 was 1431, checked in by ymipsl, 3 years ago

Local mask computation must not take into account zoom attributes.

YM

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