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

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

Updating compressed index output on using 2-level server

+) Update compressed index output with new grid distribution

Test
+) On Curie
+) test_complete:

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