source: XIOS/trunk/src/node/field.cpp @ 1119

Last change on this file since 1119 was 1119, checked in by rlacroix, 7 years ago

Add the ability to check if an output field is active at the current timestep.

The "xios_field_is_active" function can now take an optional logical argument that can be used to enable this new behavior.

  • 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: 47.0 KB
Line 
1#include "field.hpp"
2
3#include "attribute_template.hpp"
4#include "object_template.hpp"
5#include "group_template.hpp"
6
7#include "node_type.hpp"
8#include "calendar_util.hpp"
9#include "message.hpp"
10#include "xios_spl.hpp"
11#include "type.hpp"
12#include "timer.hpp"
13#include "context_client.hpp"
14#include "context_server.hpp"
15#include <set>
16#include "garbage_collector.hpp"
17#include "source_filter.hpp"
18#include "store_filter.hpp"
19#include "file_writer_filter.hpp"
20#include "pass_through_filter.hpp"
21#include "filter_expr_node.hpp"
22#include "lex_parser.hpp"
23#include "temporal_filter.hpp"
24#include "spatial_transform_filter.hpp"
25
26namespace xios{
27
28   /// ////////////////////// Définitions ////////////////////// ///
29
30   CField::CField(void)
31      : CObjectTemplate<CField>(), CFieldAttributes()
32      , grid(), file()
33      , written(false)
34      , nstep(0), nstepMax(0)
35      , hasOutputFile(false)
36      , domAxisScalarIds_(vector<StdString>(3,"")), areAllReferenceSolved(false), isReferenceSolved(false)
37      , useCompressedOutput(false)
38      , hasTimeInstant(false)
39      , hasTimeCentered(false)
40      , wasDataAlreadyReceivedFromServer(false)
41      , isEOF(false)
42   { setVirtualVariableGroup(CVariableGroup::create(getId() + "_virtual_variable_group")); }
43
44   CField::CField(const StdString& id)
45      : CObjectTemplate<CField>(id), CFieldAttributes()
46      , grid(), file()
47      , written(false)
48      , nstep(0), nstepMax(0)
49      , hasOutputFile(false)
50      , domAxisScalarIds_(vector<StdString>(3,"")), areAllReferenceSolved(false), isReferenceSolved(false)
51      , useCompressedOutput(false)
52      , hasTimeInstant(false)
53      , hasTimeCentered(false)
54      , wasDataAlreadyReceivedFromServer(false)
55      , isEOF(false)
56   { setVirtualVariableGroup(CVariableGroup::create(getId() + "_virtual_variable_group")); }
57
58   CField::~CField(void)
59   {}
60
61  //----------------------------------------------------------------
62
63   void CField::setVirtualVariableGroup(CVariableGroup* newVVariableGroup)
64   {
65      this->vVariableGroup = newVVariableGroup;
66   }
67
68   CVariableGroup* CField::getVirtualVariableGroup(void) const
69   {
70      return this->vVariableGroup;
71   }
72
73   std::vector<CVariable*> CField::getAllVariables(void) const
74   {
75      return this->vVariableGroup->getAllChildren();
76   }
77
78   void CField::solveDescInheritance(bool apply, const CAttributeMap* const parent)
79   {
80      SuperClassAttribute::setAttributes(parent, apply);
81      this->getVirtualVariableGroup()->solveDescInheritance(apply, NULL);
82   }
83
84  //----------------------------------------------------------------
85
86  bool CField::dispatchEvent(CEventServer& event)
87  {
88    if (SuperClass::dispatchEvent(event)) return true;
89    else
90    {
91      switch(event.type)
92      {
93        case EVENT_ID_UPDATE_DATA :
94          recvUpdateData(event);
95          return true;
96          break;
97
98        case EVENT_ID_READ_DATA :
99          recvReadDataRequest(event);
100          return true;
101          break;
102
103        case EVENT_ID_READ_DATA_READY :
104          recvReadDataReady(event);
105          return true;
106          break;
107
108        case EVENT_ID_ADD_VARIABLE :
109          recvAddVariable(event);
110          return true;
111          break;
112
113        case EVENT_ID_ADD_VARIABLE_GROUP :
114          recvAddVariableGroup(event);
115          return true;
116          break;
117
118        default :
119          ERROR("bool CField::dispatchEvent(CEventServer& event)", << "Unknown Event");
120          return false;
121      }
122    }
123  }
124
125  void CField::sendUpdateData(const CArray<double,1>& data)
126  {
127    CTimer::get("XIOS Send Data").resume();
128
129    CContext* context = CContext::getCurrent();
130    CContextClient* client = context->client;
131
132    CEventClient event(getType(), EVENT_ID_UPDATE_DATA);
133
134    map<int, CArray<int,1> >::iterator it;
135    list<CMessage> list_msg;
136    list<CArray<double,1> > list_data;
137
138    if (!grid->doGridHaveDataDistributed())
139    {
140       if (client->isServerLeader())
141       {
142          for (it = grid->storeIndex_toSrv.begin(); it != grid->storeIndex_toSrv.end(); it++)
143          {
144            int rank = it->first;
145            CArray<int,1>& index = it->second;
146
147            list_msg.push_back(CMessage());
148            list_data.push_back(CArray<double,1>(index.numElements()));
149
150            CArray<double,1>& data_tmp = list_data.back();
151            for (int n = 0; n < data_tmp.numElements(); n++) data_tmp(n) = data(index(n));
152
153            list_msg.back() << getId() << data_tmp;
154            event.push(rank, 1, list_msg.back());
155          }
156          client->sendEvent(event);
157       } 
158       else client->sendEvent(event);
159    }
160    else
161    {
162      for (it = grid->storeIndex_toSrv.begin(); it != grid->storeIndex_toSrv.end(); it++)
163      {
164        int rank = it->first;
165        CArray<int,1>& index = it->second;
166
167        list_msg.push_back(CMessage());
168        list_data.push_back(CArray<double,1>(index.numElements()));
169
170        CArray<double,1>& data_tmp = list_data.back();
171        for (int n = 0; n < data_tmp.numElements(); n++) data_tmp(n) = data(index(n));
172
173        list_msg.back() << getId() << data_tmp;
174        event.push(rank, grid->nbSenders[rank], list_msg.back());
175      }
176      client->sendEvent(event);
177    }
178
179    CTimer::get("XIOS Send Data").suspend();
180  }
181
182  void CField::recvUpdateData(CEventServer& event)
183  {
184    vector<int> ranks;
185    vector<CBufferIn*> buffers;
186
187    list<CEventServer::SSubEvent>::iterator it;
188    string fieldId;
189
190    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
191    {
192      int rank = it->rank;
193      CBufferIn* buffer = it->buffer;
194      *buffer >> fieldId;
195      ranks.push_back(rank);
196      buffers.push_back(buffer);
197    }
198    get(fieldId)->recvUpdateData(ranks,buffers);
199  }
200
201  void  CField::recvUpdateData(vector<int>& ranks, vector<CBufferIn*>& buffers)
202  {
203    if (data_srv.empty())
204    {
205      for (map<int, CArray<size_t, 1> >::iterator it = grid->outIndexFromClient.begin(); it != grid->outIndexFromClient.end(); ++it)
206      {
207        int rank = it->first;
208        data_srv.insert(std::make_pair(rank, CArray<double,1>(it->second.numElements())));
209        foperation_srv.insert(pair<int,boost::shared_ptr<func::CFunctor> >(rank,boost::shared_ptr<func::CFunctor>(new func::CInstant(data_srv[rank]))));
210      }
211    }
212
213    CContext* context = CContext::getCurrent();
214    const CDate& currDate = context->getCalendar()->getCurrentDate();
215    const CDate opeDate      = last_operation_srv +freq_op + freq_operation_srv - freq_op;
216    const CDate writeDate    = last_Write_srv     + freq_write_srv;
217
218    if (opeDate <= currDate)
219    {
220      for (int n = 0; n < ranks.size(); n++)
221      {
222        CArray<double,1> data_tmp;
223        *buffers[n] >> data_tmp;
224        (*foperation_srv[ranks[n]])(data_tmp);
225      }
226      last_operation_srv = currDate;
227    }
228
229    if (writeDate < (currDate + freq_operation_srv))
230    {
231      for (int n = 0; n < ranks.size(); n++)
232      {
233        this->foperation_srv[ranks[n]]->final();
234      }
235
236      last_Write_srv = writeDate;
237      writeField();
238      lastlast_Write_srv = last_Write_srv;
239    }
240  }
241
242  void CField::writeField(void)
243  {
244    if (!getRelFile()->allDomainEmpty)
245    {
246      if (grid->doGridHaveDataToWrite() || getRelFile()->type == CFile::type_attr::one_file)
247      {
248        getRelFile()->checkFile();
249        this->incrementNStep();
250        getRelFile()->getDataOutput()->writeFieldData(CField::get(this));
251      }
252    }
253  }
254
255  bool CField::sendReadDataRequest(const CDate& tsDataRequested)
256  {
257    CContext* context = CContext::getCurrent();
258    CContextClient* client = context->client;
259
260    lastDataRequestedFromServer = tsDataRequested;
261
262    if (!isEOF) // No need to send the request if we already know we are at EOF
263    {
264      CEventClient event(getType(), EVENT_ID_READ_DATA);
265      if (client->isServerLeader())
266      {
267        CMessage msg;
268        msg << getId();
269        const std::list<int>& ranks = client->getRanksServerLeader();
270        for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
271          event.push(*itRank, 1, msg);
272        client->sendEvent(event);
273      }
274      else client->sendEvent(event);
275    }
276    else
277      serverSourceFilter->signalEndOfStream(tsDataRequested);
278
279    return !isEOF;
280  }
281
282  /*!
283  Send request new data read from file if need be, that is the current data is out-of-date.
284  \return true if and only if some data was requested
285  */
286  bool CField::sendReadDataRequestIfNeeded(void)
287  {
288    const CDate& currentDate = CContext::getCurrent()->getCalendar()->getCurrentDate();
289
290    bool dataRequested = false;
291
292    while (currentDate >= lastDataRequestedFromServer)
293    {
294      info(20) << "currentDate : " << currentDate << endl ;
295      info(20) << "lastDataRequestedFromServer : " << lastDataRequestedFromServer << endl ;
296      info(20) << "file->output_freq.getValue() : " << file->output_freq.getValue() << endl ;
297      info(20) << "lastDataRequestedFromServer + file->output_freq.getValue() : " << lastDataRequestedFromServer + file->output_freq << endl ;
298
299      dataRequested |= sendReadDataRequest(lastDataRequestedFromServer + file->output_freq);
300    }
301
302    return dataRequested;
303  }
304
305  void CField::recvReadDataRequest(CEventServer& event)
306  {
307    CBufferIn* buffer = event.subEvents.begin()->buffer;
308    StdString fieldId;
309    *buffer >> fieldId;
310    get(fieldId)->recvReadDataRequest();
311  }
312
313  void CField::recvReadDataRequest(void)
314  {
315    CContext* context = CContext::getCurrent();
316    CContextClient* client = context->client;
317
318    CEventClient event(getType(), EVENT_ID_READ_DATA_READY);
319    std::list<CMessage> msgs;
320
321    bool hasData = readField();
322
323    map<int, CArray<double,1> >::iterator it;
324    if (!grid->doGridHaveDataDistributed())
325    {
326       if (client->isServerLeader())
327       {
328          if (!data_srv.empty())
329          {
330            it = data_srv.begin();
331            const std::list<int>& ranks = client->getRanksServerLeader();
332            for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
333            {
334              msgs.push_back(CMessage());
335              CMessage& msg = msgs.back();
336              msg << getId();
337              if (hasData)
338                msg << getNStep() - 1 << it->second;
339              else
340                msg << int(-1);
341              event.push(*itRank, 1, msg);
342            }
343          }
344          client->sendEvent(event);
345       } 
346       else 
347       {
348          // if (!data_srv.empty())
349          // {
350          //   it = data_srv.begin();
351          //   const std::list<int>& ranks = client->getRanksServerNotLeader();
352          //   for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
353          //   {
354          //     msgs.push_back(CMessage());
355          //     CMessage& msg = msgs.back();
356          //     msg << getId();
357          //     if (hasData)
358          //       msg << getNStep() - 1 << it->second;
359          //     else
360          //       msg << int(-1);
361          //     event.push(*itRank, 1, msg);
362          //   }
363          // }
364          client->sendEvent(event);
365       }
366    }
367    else
368    {
369      for (it = data_srv.begin(); it != data_srv.end(); it++)
370      {
371        msgs.push_back(CMessage());
372        CMessage& msg = msgs.back();
373        msg << getId();
374        if (hasData)
375          msg << getNStep() - 1 << it->second;
376        else
377          msg << int(-1);
378        event.push(it->first, grid->nbSenders[it->first], msg);
379      }
380      client->sendEvent(event);
381    }
382  }
383
384  bool CField::readField(void)
385  {
386    if (!getRelFile()->allDomainEmpty)
387    {
388      if (grid->doGridHaveDataToWrite() || getRelFile()->type == CFile::type_attr::one_file)
389      {
390        if (data_srv.empty())
391        {
392          for (map<int, CArray<size_t, 1> >::iterator it = grid->outIndexFromClient.begin(); it != grid->outIndexFromClient.end(); ++it)
393            data_srv.insert(std::make_pair(it->first, CArray<double,1>(it->second.numElements())));
394        }
395
396        getRelFile()->checkFile();
397        if (!nstepMax)
398        {
399          nstepMax = getRelFile()->getDataInput()->getFieldNbRecords(CField::get(this));
400        }
401
402        this->incrementNStep();
403
404        if (getNStep() > nstepMax && (getRelFile()->cyclic.isEmpty() || !getRelFile()->cyclic) )
405          return false;
406
407        getRelFile()->getDataInput()->readFieldData(CField::get(this));
408      }
409    }
410
411    return true;
412  }
413
414  void CField::recvReadDataReady(CEventServer& event)
415  {
416    string fieldId;
417    vector<int> ranks;
418    vector<CBufferIn*> buffers;
419
420    list<CEventServer::SSubEvent>::iterator it;
421    for (it = event.subEvents.begin(); it != event.subEvents.end(); ++it)
422    {
423      ranks.push_back(it->rank);
424      CBufferIn* buffer = it->buffer;
425      *buffer >> fieldId;
426      buffers.push_back(buffer);
427    }
428    get(fieldId)->recvReadDataReady(ranks, buffers);
429  }
430
431  void CField::recvReadDataReady(vector<int> ranks, vector<CBufferIn*> buffers)
432  {
433    CContext* context = CContext::getCurrent();
434    int record;
435    std::map<int, CArray<double,1> > data;
436
437    for (int i = 0; i < ranks.size(); i++)
438    {
439      int rank = ranks[i];
440      *buffers[i] >> record;
441      isEOF = (record == int(-1));
442
443      if (!isEOF)
444        *buffers[i] >> data[rank];
445      else
446        break;
447    }
448
449    if (wasDataAlreadyReceivedFromServer)
450      lastDataReceivedFromServer = lastDataReceivedFromServer + file->output_freq;
451    else
452    {
453      lastDataReceivedFromServer = context->getCalendar()->getInitDate();
454      wasDataAlreadyReceivedFromServer = true;
455    }
456
457    if (isEOF)
458      serverSourceFilter->signalEndOfStream(lastDataReceivedFromServer);
459    else
460      serverSourceFilter->streamDataFromServer(lastDataReceivedFromServer, data);
461  }
462
463   //----------------------------------------------------------------
464
465   void CField::setRelFile(CFile* _file)
466   {
467      this->file = _file;
468      hasOutputFile = true;
469   }
470
471   //----------------------------------------------------------------
472
473   StdString CField::GetName(void)    { return StdString("field"); }
474   StdString CField::GetDefName(void) { return CField::GetName(); }
475   ENodeType CField::GetType(void)    { return eField; }
476
477   //----------------------------------------------------------------
478
479   CGrid* CField::getRelGrid(void) const
480   {
481      return this->grid;
482   }
483
484   //----------------------------------------------------------------
485
486   CFile* CField::getRelFile(void) const
487   {
488      return this->file;
489   }
490
491   int CField::getNStep(void) const
492   {
493      return this->nstep;
494   }
495
496   func::CFunctor::ETimeType CField::getOperationTimeType() const
497   {
498     return operationTimeType;
499   }
500
501   //----------------------------------------------------------------
502
503   void CField::incrementNStep(void)
504   {
505      this->nstep++;
506   }
507
508   void CField::resetNStep(int nstep /*= 0*/)
509   {
510      this->nstep = nstep;
511   }
512
513   void CField::resetNStepMax(void)
514   {
515      this->nstepMax = 0;
516   }
517
518   //----------------------------------------------------------------
519
520   bool CField::isActive(bool atCurrentTimestep /*= false*/) const
521   {
522      if (atCurrentTimestep && clientSourceFilter)
523        return clientSourceFilter->isDataExpected(CContext::getCurrent()->getCalendar()->getCurrentDate());
524      else
525        return (instantDataFilter != NULL);
526   }
527
528   //----------------------------------------------------------------
529
530   bool CField::wasWritten() const
531   {
532     return written;
533   }
534
535   void CField::setWritten()
536   {
537     written = true;
538   }
539
540   //----------------------------------------------------------------
541
542   bool CField::getUseCompressedOutput() const
543   {
544     return useCompressedOutput;
545   }
546
547   void CField::setUseCompressedOutput()
548   {
549     useCompressedOutput = true;
550   }
551
552   //----------------------------------------------------------------
553
554   boost::shared_ptr<COutputPin> CField::getInstantDataFilter()
555   {
556     return instantDataFilter;
557   }
558
559   //----------------------------------------------------------------
560
561   void CField::solveOnlyReferenceEnabledField(bool doSending2Server)
562   {
563     CContext* context = CContext::getCurrent();
564     if (!isReferenceSolved)
565     {
566        isReferenceSolved = true;
567
568        if (context->hasClient)
569        {
570          solveRefInheritance(true);
571          if (hasDirectFieldReference()) getDirectFieldReference()->solveOnlyReferenceEnabledField(false);
572        }
573        else if (context->hasServer)
574          solveServerOperation();
575
576        solveGridReference();
577
578       if (context->hasClient)
579       {
580         solveGenerateGrid();
581         buildGridTransformationGraph();
582       }
583     }
584   }
585
586   /*!
587     Build up graph of grids which plays role of destination and source in grid transformation
588     This function should be called before \func solveGridReference()
589   */
590   void CField::buildGridTransformationGraph()
591   {
592     CContext* context = CContext::getCurrent();
593     if (context->hasClient)
594     {
595       if (grid && !grid->isTransformed() && hasDirectFieldReference() && grid != getDirectFieldReference()->grid)
596       {
597         grid->addTransGridSource(getDirectFieldReference()->grid);
598       }
599     }
600   }
601
602   /*!
603     Generate a new grid destination if there are more than one grid source pointing to a same grid destination
604   */
605   void CField::generateNewTransformationGridDest()
606   {
607     CContext* context = CContext::getCurrent();
608     if (context->hasClient)
609     {
610       std::map<CGrid*,std::pair<bool,StdString> >& gridSrcMap = grid->getTransGridSource();
611       if (1 < gridSrcMap.size())
612       {
613         // Search for grid source
614         CGrid* gridSrc = grid;
615         CField* currField = this;
616         std::vector<CField*> hieraField;
617
618         while (currField->hasDirectFieldReference() && (gridSrc == grid))
619         {
620           hieraField.push_back(currField);
621           CField* tmp = currField->getDirectFieldReference();
622           currField = tmp;
623           gridSrc = currField->grid;
624         }
625
626         if (gridSrcMap.end() != gridSrcMap.find(gridSrc))
627         {
628           CGrid* gridTmp;
629           std::pair<bool,StdString> newGridDest = gridSrcMap[gridSrc];
630           if (newGridDest.first)
631           {
632             StdString newIdGridDest = newGridDest.second;
633             if (!CGrid::has(newIdGridDest))
634             {
635                ERROR("CGrid* CGrid::generateNewTransformationGridDest()",
636                  << " Something wrong happened! Grid whose id " << newIdGridDest
637                  << "should exist ");
638             }
639             gridTmp = CGrid::get(newIdGridDest);
640           }
641           else
642           {
643             StdString newIdGridDest = CGrid::generateId(gridSrc, grid);
644             gridTmp = CGrid::cloneGrid(newIdGridDest, grid);
645
646             (gridSrcMap[gridSrc]).first = true;
647             (gridSrcMap[gridSrc]).second = newIdGridDest;
648           }
649
650           // Update all descendants
651           for (std::vector<CField*>::iterator it = hieraField.begin(); it != hieraField.end(); ++it)
652           {
653             (*it)->grid = gridTmp;
654             (*it)->updateRef((*it)->grid);
655           }
656         }
657       }
658     }
659   }
660
661   void CField::updateRef(CGrid* grid)
662   {
663     if (!grid_ref.isEmpty()) grid_ref.setValue(grid->getId());
664     else
665     {
666       std::vector<CAxis*> axisTmp = grid->getAxis();
667       std::vector<CDomain*> domainTmp = grid->getDomains();
668       if ((1<axisTmp.size()) || (1<domainTmp.size()))
669         ERROR("void CField::updateRef(CGrid* grid)",
670           << "More than one domain or axis is available for domain_ref/axis_ref of field " << this->getId());
671
672       if ((!domain_ref.isEmpty()) && (domainTmp.empty()))
673         ERROR("void CField::updateRef(CGrid* grid)",
674           << "Incoherent between available domain and domain_ref of field " << this->getId());
675       if ((!axis_ref.isEmpty()) && (axisTmp.empty()))
676         ERROR("void CField::updateRef(CGrid* grid)",
677           << "Incoherent between available axis and axis_ref of field " << this->getId());
678
679       if (!domain_ref.isEmpty()) domain_ref.setValue(domainTmp[0]->getId());
680       if (!axis_ref.isEmpty()) axis_ref.setValue(axisTmp[0]->getId());
681     }
682   }
683
684   void CField::solveAllReferenceEnabledField(bool doSending2Server)
685   {
686     CContext* context = CContext::getCurrent();
687     solveOnlyReferenceEnabledField(doSending2Server);
688
689     if (!areAllReferenceSolved)
690     {
691        areAllReferenceSolved = true;
692
693        if (context->hasClient)
694        {
695          solveRefInheritance(true);
696          if (hasDirectFieldReference()) getDirectFieldReference()->solveAllReferenceEnabledField(false);
697        }
698        else if (context->hasServer)
699          solveServerOperation();
700
701        solveGridReference();
702     }
703
704     solveGridDomainAxisRef(doSending2Server);
705
706     if (context->hasClient)
707     {
708       solveTransformedGrid();
709     }
710
711     solveCheckMaskIndex(doSending2Server);
712   }
713
714   std::map<int, StdSize> CField::getGridAttributesBufferSize()
715   {
716     return grid->getAttributesBufferSize();
717   }
718
719   std::map<int, StdSize> CField::getGridDataBufferSize()
720   {
721     return grid->getDataBufferSize(getId());
722   }
723
724   //----------------------------------------------------------------
725
726   void CField::solveServerOperation(void)
727   {
728      CContext* context = CContext::getCurrent();
729
730      if (!context->hasServer || !hasOutputFile) return;
731
732      if (freq_op.isEmpty())
733        freq_op.setValue(TimeStep);
734
735      if (freq_offset.isEmpty())
736        freq_offset.setValue(NoneDu);
737
738      freq_operation_srv = file->output_freq.getValue();
739      freq_write_srv     = file->output_freq.getValue();
740
741      lastlast_Write_srv = context->getCalendar()->getInitDate();
742      last_Write_srv     = context->getCalendar()->getInitDate();
743      last_operation_srv = context->getCalendar()->getInitDate();
744
745      const CDuration toffset = freq_operation_srv - freq_offset.getValue() - context->getCalendar()->getTimeStep();
746      last_operation_srv     = last_operation_srv - toffset;
747
748      if (operation.isEmpty())
749        ERROR("void CField::solveServerOperation(void)",
750              << "An operation must be defined for field \"" << getId() << "\".");
751
752      boost::shared_ptr<func::CFunctor> functor;
753      CArray<double, 1> dummyData;
754
755#define DECLARE_FUNCTOR(MType, mtype) \
756      if (operation.getValue().compare(#mtype) == 0) \
757      { \
758        functor.reset(new func::C##MType(dummyData)); \
759      }
760
761#include "functor_type.conf"
762
763      if (!functor)
764        ERROR("void CField::solveServerOperation(void)",
765              << "\"" << operation << "\" is not a valid operation.");
766
767      operationTimeType = functor->timeType();
768   }
769
770   //----------------------------------------------------------------
771
772   /*!
773    * Constructs the graph filter for the field, enabling or not the data output.
774    * This method should not be called more than once with enableOutput equal to true.
775    *
776    * \param gc the garbage collector to use when building the filter graph
777    * \param enableOutput must be true when the field data is to be
778    *                     read by the client or/and written to a file
779    */
780   void CField::buildFilterGraph(CGarbageCollector& gc, bool enableOutput)
781   {
782     if (!areAllReferenceSolved) solveAllReferenceEnabledField(false);
783
784     // Start by building a filter which can provide the field's instant data
785     if (!instantDataFilter)
786     {
787       // Check if we have an expression to parse
788       if (hasExpression())
789       {
790         boost::scoped_ptr<IFilterExprNode> expr(parseExpr(getExpression() + '\0'));
791         boost::shared_ptr<COutputPin> filter = expr->reduce(gc, *this);
792
793         // Check if a spatial transformation is needed
794         if (!field_ref.isEmpty())
795         {
796           CGrid* gridRef = CField::get(field_ref)->grid;
797
798           if (grid && grid != gridRef && grid->hasTransform())
799           {
800             bool hasMissingValue = (!detect_missing_value.isEmpty() && !default_value.isEmpty() && detect_missing_value == true);
801             double defaultValue  = hasMissingValue ? default_value : (!default_value.isEmpty() ? default_value : 0.0);
802             std::pair<boost::shared_ptr<CFilter>, boost::shared_ptr<CFilter> > filters = CSpatialTransformFilter::buildFilterGraph(gc, gridRef, grid, hasMissingValue, defaultValue);
803
804             filter->connectOutput(filters.first, 0);
805             filter = filters.second;
806           }
807         }
808
809         instantDataFilter = filter;
810       }
811       // Check if we have a reference on another field
812       else if (!field_ref.isEmpty())
813         instantDataFilter = getFieldReference(gc);
814       // Check if the data is to be read from a file
815       else if (file && !file->mode.isEmpty() && file->mode == CFile::mode_attr::read)
816         instantDataFilter = serverSourceFilter = boost::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid,
817                                                                                                     freq_offset.isEmpty() ? NoneDu : freq_offset,
818                                                                                                     true));
819       else // The data might be passed from the model
820       {
821          bool ignoreMissingValue = (!detect_missing_value.isEmpty() && !default_value.isEmpty() && detect_missing_value == true);
822          double defaultValue  = ignoreMissingValue ? default_value : (!default_value.isEmpty() ? default_value : 0.0);
823          instantDataFilter = clientSourceFilter = boost::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, NoneDu, false,
824                                                                                                      ignoreMissingValue, defaultValue));       }
825     }
826
827     // If the field data is to be read by the client or/and written to a file
828     if (enableOutput && !storeFilter && !fileWriterFilter)
829     {
830       if (!read_access.isEmpty() && read_access)
831       {
832         storeFilter = boost::shared_ptr<CStoreFilter>(new CStoreFilter(gc, CContext::getCurrent(), grid));
833         instantDataFilter->connectOutput(storeFilter, 0);
834       }
835
836       if (file && (file->mode.isEmpty() || file->mode == CFile::mode_attr::write))
837       {
838         fileWriterFilter = boost::shared_ptr<CFileWriterFilter>(new CFileWriterFilter(gc, this));
839         getTemporalDataFilter(gc, file->output_freq)->connectOutput(fileWriterFilter, 0);
840       }
841     }
842   }
843
844   /*!
845    * Returns the filter needed to handle the field reference.
846    * This method should only be called when building the filter graph corresponding to the field.
847    *
848    * \param gc the garbage collector to use
849    * \return the output pin corresponding to the field reference
850    */
851   boost::shared_ptr<COutputPin> CField::getFieldReference(CGarbageCollector& gc)
852   {
853     if (instantDataFilter || field_ref.isEmpty())
854       ERROR("COutputPin* CField::getFieldReference(CGarbageCollector& gc)",
855             "Impossible to get the field reference for a field which has already been parsed or which does not have a field_ref.");
856
857     CField* fieldRef = CField::get(field_ref);
858     fieldRef->buildFilterGraph(gc, false);
859
860     std::pair<boost::shared_ptr<CFilter>, boost::shared_ptr<CFilter> > filters;
861     // Check if a spatial transformation is needed
862     if (grid && grid != fieldRef->grid && grid->hasTransform())
863     {       
864       bool hasMissingValue = (!detect_missing_value.isEmpty() && !default_value.isEmpty() && detect_missing_value == true);
865       double defaultValue  = hasMissingValue ? default_value : (!default_value.isEmpty() ? default_value : 0.0);                               
866       filters = CSpatialTransformFilter::buildFilterGraph(gc, fieldRef->grid, grid, hasMissingValue, defaultValue);
867     }
868     else
869       filters.first = filters.second = boost::shared_ptr<CFilter>(new CPassThroughFilter(gc));
870
871     fieldRef->getInstantDataFilter()->connectOutput(filters.first, 0);
872
873     return filters.second;
874   }
875
876   /*!
877    * Returns the filter needed to handle a self reference in the field's expression.
878    * If the needed filter does not exist, it is created, otherwise it is reused.
879    * This method should only be called when building the filter graph corresponding
880    * to the field's expression.
881    *
882    * \param gc the garbage collector to use
883    * \return the output pin corresponding to a self reference
884    */
885   boost::shared_ptr<COutputPin> CField::getSelfReference(CGarbageCollector& gc)
886   {
887     if (instantDataFilter || !hasExpression())
888       ERROR("COutputPin* CField::getSelfReference(CGarbageCollector& gc)",
889             "Impossible to add a self reference to a field which has already been parsed or which does not have an expression.");
890
891     if (!selfReferenceFilter)
892     {
893       if (file && !file->mode.isEmpty() && file->mode == CFile::mode_attr::read)
894       {
895         if (!serverSourceFilter)
896           serverSourceFilter = boost::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid,
897                                                                                   freq_offset.isEmpty() ? NoneDu : freq_offset,
898                                                                                   true));
899
900         selfReferenceFilter = serverSourceFilter;
901       }
902       else if (!field_ref.isEmpty())
903       {
904         CField* fieldRef = CField::get(field_ref);
905         fieldRef->buildFilterGraph(gc, false); 
906         selfReferenceFilter = fieldRef->getInstantDataFilter();
907       }
908       else
909       {
910         if (!clientSourceFilter)
911         {
912           bool ignoreMissingValue = (!detect_missing_value.isEmpty() && !default_value.isEmpty() && detect_missing_value == true);
913           double defaultValue  = ignoreMissingValue ? default_value : (!default_value.isEmpty() ? default_value : 0.0); 
914           clientSourceFilter = boost::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, NoneDu, false,
915                                                                                   ignoreMissingValue, defaultValue));
916         }
917
918         selfReferenceFilter = clientSourceFilter;
919       }
920     }
921
922     return selfReferenceFilter;
923   }
924
925   /*!
926    * Returns the temporal filter corresponding to the field's temporal operation
927    * for the specified operation frequency. The filter is created if it does not
928    * exist, otherwise it is reused.
929    *
930    * \param gc the garbage collector to use
931    * \param outFreq the operation frequency, i.e. the frequency at which the output data will be computed
932    * \return the output pin corresponding to the requested temporal filter
933    */
934   boost::shared_ptr<COutputPin> CField::getTemporalDataFilter(CGarbageCollector& gc, CDuration outFreq)
935   {
936     std::map<CDuration, boost::shared_ptr<COutputPin> >::iterator it = temporalDataFilters.find(outFreq);
937
938     if (it == temporalDataFilters.end())
939     {
940       if (operation.isEmpty())
941         ERROR("void CField::getTemporalDataFilter(CGarbageCollector& gc, CDuration outFreq)",
942               << "An operation must be defined for field \"" << getId() << "\".");
943
944       if (freq_op.isEmpty())
945         freq_op.setValue(TimeStep);
946       if (freq_offset.isEmpty())
947         freq_offset.setValue(NoneDu);
948
949       const bool ignoreMissingValue = (!detect_missing_value.isEmpty() && !default_value.isEmpty() && detect_missing_value == true);
950       
951       boost::shared_ptr<CTemporalFilter> temporalFilter(new CTemporalFilter(gc, operation,
952                                                                             CContext::getCurrent()->getCalendar()->getInitDate(),
953                                                                             freq_op, freq_offset, outFreq,
954                                                                             ignoreMissingValue, ignoreMissingValue ? default_value : 0.0));
955       instantDataFilter->connectOutput(temporalFilter, 0);
956
957       it = temporalDataFilters.insert(std::make_pair(outFreq, temporalFilter)).first;
958     }
959
960     return it->second;
961   }
962
963  /*!
964    * Returns the temporal filter corresponding to the field's temporal operation
965    * for the specified operation frequency.
966    *
967    * \param gc the garbage collector to use
968    * \param outFreq the operation frequency, i.e. the frequency at which the output data will be computed
969    * \return the output pin corresponding to the requested temporal filter
970    */
971   
972   boost::shared_ptr<COutputPin> CField::getSelfTemporalDataFilter(CGarbageCollector& gc, CDuration outFreq)
973   {
974     if (instantDataFilter || !hasExpression())
975       ERROR("COutputPin* CField::getSelfTemporalDataFilter(CGarbageCollector& gc)",
976             "Impossible to add a self reference to a field which has already been parsed or which does not have an expression.");
977
978     if (!selfReferenceFilter) getSelfReference(gc) ;
979
980     if (serverSourceFilter || clientSourceFilter)
981     {
982       if (operation.isEmpty())
983         ERROR("void CField::getSelfTemporalDataFilter(CGarbageCollector& gc, CDuration outFreq)",
984               << "An operation must be defined for field \"" << getId() << "\".");
985
986       if (freq_op.isEmpty()) freq_op.setValue(TimeStep);
987       if (freq_offset.isEmpty()) freq_offset.setValue(NoneDu);
988
989       const bool ignoreMissingValue = (!detect_missing_value.isEmpty() && !default_value.isEmpty() && detect_missing_value == true);
990
991       boost::shared_ptr<CTemporalFilter> temporalFilter(new CTemporalFilter(gc, operation,
992                                                                             CContext::getCurrent()->getCalendar()->getInitDate(),
993                                                                             freq_op, freq_offset, outFreq,
994                                                                             ignoreMissingValue, ignoreMissingValue ? default_value : 0.0));
995       selfReferenceFilter->connectOutput(temporalFilter, 0);
996       return temporalFilter ;
997     }
998     else if (!field_ref.isEmpty())
999     {
1000       CField* fieldRef = CField::get(field_ref);
1001       fieldRef->buildFilterGraph(gc, false); 
1002       return fieldRef->getTemporalDataFilter(gc, outFreq) ;
1003     }
1004  }
1005
1006   //----------------------------------------------------------------
1007/*
1008   void CField::fromBinary(StdIStream& is)
1009   {
1010      SuperClass::fromBinary(is);
1011#define CLEAR_ATT(name_)\
1012      SuperClassAttribute::operator[](#name_)->reset()
1013
1014         CLEAR_ATT(domain_ref);
1015         CLEAR_ATT(axis_ref);
1016#undef CLEAR_ATT
1017
1018   }
1019*/
1020   //----------------------------------------------------------------
1021
1022   void CField::solveGridReference(void)
1023   {
1024      if (grid_ref.isEmpty() && domain_ref.isEmpty() && axis_ref.isEmpty() && scalar_ref.isEmpty())
1025      {
1026        ERROR("CField::solveGridReference(void)",
1027              << "A grid must be defined for field '" << getFieldOutputName() << "' .");
1028      }
1029      else if (!grid_ref.isEmpty() && (!domain_ref.isEmpty() || !axis_ref.isEmpty() || !scalar_ref.isEmpty()))
1030      {
1031        ERROR("CField::solveGridReference(void)",
1032              << "Field '" << getFieldOutputName() << "' has both a grid and a domain/axis/scalar." << std::endl
1033              << "Please define either 'grid_ref' or 'domain_ref'/'axis_ref'/'scalar_ref'.");
1034      }
1035
1036      if (grid_ref.isEmpty())
1037      {
1038        std::vector<CDomain*> vecDom;
1039        std::vector<CAxis*> vecAxis;
1040        std::vector<CScalar*> vecScalar;
1041        std::vector<int> axisDomainOrderTmp;
1042       
1043        if (!domain_ref.isEmpty())
1044        {
1045          StdString tmp = domain_ref.getValue();
1046          if (CDomain::has(domain_ref))
1047          {
1048            vecDom.push_back(CDomain::get(domain_ref));
1049            axisDomainOrderTmp.push_back(2);
1050          }
1051          else
1052            ERROR("CField::solveGridReference(void)",
1053                  << "Invalid reference to domain '" << domain_ref.getValue() << "'.");
1054        }
1055
1056        if (!axis_ref.isEmpty())
1057        {
1058          if (CAxis::has(axis_ref))
1059          {
1060            vecAxis.push_back(CAxis::get(axis_ref));
1061            axisDomainOrderTmp.push_back(1);
1062          }
1063          else
1064            ERROR("CField::solveGridReference(void)",
1065                  << "Invalid reference to axis '" << axis_ref.getValue() << "'.");
1066        }
1067
1068        if (!scalar_ref.isEmpty())
1069        {
1070          if (CScalar::has(scalar_ref))
1071          {
1072            vecScalar.push_back(CScalar::get(scalar_ref));
1073            axisDomainOrderTmp.push_back(0);
1074          }
1075          else
1076            ERROR("CField::solveGridReference(void)",
1077                  << "Invalid reference to scalar '" << scalar_ref.getValue() << "'.");
1078        }
1079       
1080        CArray<int,1> axisDomainOrder(axisDomainOrderTmp.size());
1081        for (int idx = 0; idx < axisDomainOrderTmp.size(); ++idx)
1082        {
1083          axisDomainOrder(idx) = axisDomainOrderTmp[idx];
1084        }
1085
1086        // Warning: the gridId shouldn't be set as the grid_ref since it could be inherited
1087        StdString gridId = CGrid::generateId(vecDom, vecAxis, vecScalar,axisDomainOrder);
1088        if (CGrid::has(gridId))
1089          this->grid = CGrid::get(gridId);
1090        else
1091          this->grid = CGrid::createGrid(gridId, vecDom, vecAxis, vecScalar,axisDomainOrder);
1092      }
1093      else
1094      {
1095        if (CGrid::has(grid_ref))
1096          this->grid = CGrid::get(grid_ref);
1097        else
1098          ERROR("CField::solveGridReference(void)",
1099                << "Invalid reference to grid '" << grid_ref.getValue() << "'.");
1100      }
1101   }
1102
1103   void CField::solveGridDomainAxisRef(bool checkAtt)
1104   {
1105     grid->solveDomainAxisRef(checkAtt);
1106   }
1107
1108   void CField::solveCheckMaskIndex(bool doSendingIndex)
1109   {
1110     grid->checkMaskIndex(doSendingIndex);
1111   }
1112
1113   void CField::solveTransformedGrid()
1114   {
1115     if (grid && !grid->isTransformed() && hasDirectFieldReference() && grid != getDirectFieldReference()->grid)
1116     {
1117       std::vector<CGrid*> grids;
1118       // Source grid
1119       grids.push_back(getDirectFieldReference()->grid);
1120       // Intermediate grids
1121       if (!grid_path.isEmpty())
1122       {
1123         std::string gridId;
1124         size_t start = 0, end;
1125
1126         do
1127         {
1128           end = grid_path.getValue().find(',', start);
1129           if (end != std::string::npos)
1130           {
1131             gridId = grid_path.getValue().substr(start, end - start);
1132             start = end + 1;
1133           }
1134           else
1135             gridId = grid_path.getValue().substr(start);
1136
1137           if (!CGrid::has(gridId))
1138             ERROR("void CField::solveTransformedGrid()",
1139                   << "Invalid grid_path, the grid '" << gridId << "' does not exist.");
1140
1141           grids.push_back(CGrid::get(gridId));
1142         }
1143         while (end != std::string::npos);
1144       }
1145       // Destination grid
1146       grids.push_back(grid);
1147
1148       for (size_t i = 0, count = grids.size() - 1; i < count; ++i)
1149       {
1150         CGrid *gridSrc  = grids[i];
1151         CGrid *gridDest = grids[i + 1];
1152         if (!gridDest->isTransformed())
1153           gridDest->transformGrid(gridSrc);
1154       }
1155     }
1156     else if (grid && grid->hasTransform() && !grid->isTransformed())
1157     {
1158       // Temporarily deactivate the self-transformation of grid
1159       //grid->transformGrid(grid);
1160     }
1161   }
1162
1163   void CField::solveGenerateGrid()
1164   {
1165     if (grid && !grid->isTransformed() && hasDirectFieldReference() && grid != getDirectFieldReference()->grid)
1166       grid->completeGrid(getDirectFieldReference()->grid);
1167     else
1168       grid->completeGrid();
1169   }
1170
1171   void CField::solveGridDomainAxisBaseRef()
1172   {
1173     grid->solveDomainAxisRef(false);
1174     grid->solveDomainAxisBaseRef();
1175   }
1176
1177   ///-------------------------------------------------------------------
1178
1179   template <>
1180   void CGroupTemplate<CField, CFieldGroup, CFieldAttributes>::solveRefInheritance(void)
1181   {
1182      if (this->group_ref.isEmpty()) return;
1183      StdString gref = this->group_ref.getValue();
1184
1185      if (!CFieldGroup::has(gref))
1186         ERROR("CGroupTemplate<CField, CFieldGroup, CFieldAttributes>::solveRefInheritance(void)",
1187               << "[ gref = " << gref << "]"
1188               << " invalid group name !");
1189
1190      CFieldGroup* group = CFieldGroup::get(gref);
1191      CFieldGroup* owner = CFieldGroup::get(boost::polymorphic_downcast<CFieldGroup*>(this));
1192
1193      std::vector<CField*> allChildren  = group->getAllChildren();
1194      std::vector<CField*>::iterator it = allChildren.begin(), end = allChildren.end();
1195
1196      for (; it != end; it++)
1197      {
1198         CField* child = *it;
1199         if (child->hasId()) owner->createChild()->field_ref.setValue(child->getId());
1200
1201      }
1202   }
1203
1204   void CField::scaleFactorAddOffset(double scaleFactor, double addOffset)
1205   {
1206     map<int, CArray<double,1> >::iterator it;
1207     for (it = data_srv.begin(); it != data_srv.end(); it++) it->second = (it->second - addOffset) / scaleFactor;
1208   }
1209
1210   void CField::invertScaleFactorAddOffset(double scaleFactor, double addOffset)
1211   {
1212     map<int, CArray<double,1> >::iterator it;
1213     for (it = data_srv.begin(); it != data_srv.end(); it++) it->second = it->second * scaleFactor + addOffset;
1214   }
1215
1216   void CField::outputField(CArray<double,3>& fieldOut)
1217   {
1218      map<int, CArray<double,1> >::iterator it;
1219      for (it = data_srv.begin(); it != data_srv.end(); it++)
1220      {
1221        grid->outputField(it->first, it->second, fieldOut.dataFirst());
1222      }
1223   }
1224
1225   void CField::outputField(CArray<double,2>& fieldOut)
1226   {
1227      map<int, CArray<double,1> >::iterator it;
1228      for(it=data_srv.begin();it!=data_srv.end();it++)
1229      {
1230         grid->outputField(it->first, it->second, fieldOut.dataFirst());
1231      }
1232   }
1233
1234   void CField::outputField(CArray<double,1>& fieldOut)
1235   {
1236      map<int, CArray<double,1> >::iterator it;
1237
1238      for (it = data_srv.begin(); it != data_srv.end(); it++)
1239      {
1240         grid->outputField(it->first, it->second, fieldOut.dataFirst());
1241      }
1242   }
1243
1244   void CField::inputField(CArray<double,3>& fieldOut)
1245   {
1246      map<int, CArray<double,1> >::iterator it;
1247      for (it = data_srv.begin(); it != data_srv.end(); it++)
1248      {
1249        grid->inputField(it->first, fieldOut.dataFirst(), it->second);
1250      }
1251   }
1252
1253   void CField::inputField(CArray<double,2>& fieldOut)
1254   {
1255      map<int, CArray<double,1> >::iterator it;
1256      for(it = data_srv.begin(); it != data_srv.end(); it++)
1257      {
1258         grid->inputField(it->first, fieldOut.dataFirst(), it->second);
1259      }
1260   }
1261
1262   void CField::inputField(CArray<double,1>& fieldOut)
1263   {
1264      map<int, CArray<double,1> >::iterator it;
1265      for (it = data_srv.begin(); it != data_srv.end(); it++)
1266      {
1267         grid->inputField(it->first, fieldOut.dataFirst(), it->second);
1268      }
1269   }
1270
1271   void CField::outputCompressedField(CArray<double,1>& fieldOut)
1272   {
1273      map<int, CArray<double,1> >::iterator it;
1274
1275      for (it = data_srv.begin(); it != data_srv.end(); it++)
1276      {
1277         grid->outputCompressedField(it->first, it->second, fieldOut.dataFirst());
1278      }
1279   }
1280
1281   ///-------------------------------------------------------------------
1282
1283   void CField::parse(xml::CXMLNode& node)
1284   {
1285      SuperClass::parse(node);
1286      if (!node.getContent(this->content))
1287      {
1288        if (node.goToChildElement())
1289        {
1290          do
1291          {
1292            if (node.getElementName() == "variable" || node.getElementName() == "variable_group") this->getVirtualVariableGroup()->parseChild(node);
1293          } while (node.goToNextElement());
1294          node.goToParentElement();
1295        }
1296      }
1297    }
1298
1299   /*!
1300     This function retrieves Id of corresponding domain_ref and axis_ref (if any)
1301   of a field. In some cases, only domain exists but axis doesn't
1302   \return pair of Domain and Axis id
1303   */
1304   const std::vector<StdString>& CField::getRefDomainAxisIds()
1305   {
1306     CGrid* cgPtr = getRelGrid();
1307     if (NULL != cgPtr)
1308     {
1309       std::vector<StdString>::iterator it;
1310       if (!domain_ref.isEmpty())
1311       {
1312         std::vector<StdString> domainList = cgPtr->getDomainList();
1313         it = std::find(domainList.begin(), domainList.end(), domain_ref.getValue());
1314         if (domainList.end() != it) domAxisScalarIds_[0] = *it;
1315       }
1316
1317       if (!axis_ref.isEmpty())
1318       {
1319         std::vector<StdString> axisList = cgPtr->getAxisList();
1320         it = std::find(axisList.begin(), axisList.end(), axis_ref.getValue());
1321         if (axisList.end() != it) domAxisScalarIds_[1] = *it;
1322       }
1323
1324       if (!scalar_ref.isEmpty())
1325       {
1326         std::vector<StdString> scalarList = cgPtr->getScalarList();
1327         it = std::find(scalarList.begin(), scalarList.end(), scalar_ref.getValue());
1328         if (scalarList.end() != it) domAxisScalarIds_[2] = *it;
1329       }
1330     }
1331     return (domAxisScalarIds_);
1332   }
1333
1334   CVariable* CField::addVariable(const string& id)
1335   {
1336     return vVariableGroup->createChild(id);
1337   }
1338
1339   CVariableGroup* CField::addVariableGroup(const string& id)
1340   {
1341     return vVariableGroup->createChildGroup(id);
1342   }
1343
1344   void CField::sendAddAllVariables()
1345   {
1346     std::vector<CVariable*> allVar = getAllVariables();
1347     std::vector<CVariable*>::const_iterator it = allVar.begin();
1348     std::vector<CVariable*>::const_iterator itE = allVar.end();
1349
1350     for (; it != itE; ++it)
1351     {
1352       this->sendAddVariable((*it)->getId());
1353       (*it)->sendAllAttributesToServer();
1354       (*it)->sendValue();
1355     }
1356   }
1357
1358   void CField::sendAddVariable(const string& id)
1359   {
1360    CContext* context = CContext::getCurrent();
1361
1362    if (!context->hasServer)
1363    {
1364       CContextClient* client = context->client;
1365
1366       CEventClient event(this->getType(),EVENT_ID_ADD_VARIABLE);
1367       if (client->isServerLeader())
1368       {
1369         CMessage msg;
1370         msg << this->getId();
1371         msg << id;
1372         const std::list<int>& ranks = client->getRanksServerLeader();
1373         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1374           event.push(*itRank,1,msg);
1375         client->sendEvent(event);
1376       }
1377       else client->sendEvent(event);
1378    }
1379   }
1380
1381   void CField::sendAddVariableGroup(const string& id)
1382   {
1383    CContext* context = CContext::getCurrent();
1384    if (!context->hasServer)
1385    {
1386       CContextClient* client = context->client;
1387
1388       CEventClient event(this->getType(),EVENT_ID_ADD_VARIABLE_GROUP);
1389       if (client->isServerLeader())
1390       {
1391         CMessage msg;
1392         msg << this->getId();
1393         msg << id;
1394         const std::list<int>& ranks = client->getRanksServerLeader();
1395         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1396           event.push(*itRank,1,msg);
1397         client->sendEvent(event);
1398       }
1399       else client->sendEvent(event);
1400    }
1401   }
1402
1403   void CField::recvAddVariable(CEventServer& event)
1404   {
1405
1406      CBufferIn* buffer = event.subEvents.begin()->buffer;
1407      string id;
1408      *buffer >> id;
1409      get(id)->recvAddVariable(*buffer);
1410   }
1411
1412   void CField::recvAddVariable(CBufferIn& buffer)
1413   {
1414      string id;
1415      buffer >> id;
1416      addVariable(id);
1417   }
1418
1419   void CField::recvAddVariableGroup(CEventServer& event)
1420   {
1421
1422      CBufferIn* buffer = event.subEvents.begin()->buffer;
1423      string id;
1424      *buffer >> id;
1425      get(id)->recvAddVariableGroup(*buffer);
1426   }
1427
1428   void CField::recvAddVariableGroup(CBufferIn& buffer)
1429   {
1430      string id;
1431      buffer >> id;
1432      addVariableGroup(id);
1433   }
1434
1435   /*!
1436    * Returns string arithmetic expression associated to the field.
1437    * \return if content is defined return content string, otherwise, if "expr" attribute is defined, return expr string.
1438    */
1439   const string& CField::getExpression(void)
1440   {
1441     if (!expr.isEmpty() && content.empty())
1442     {
1443       content = expr;
1444       expr.reset();
1445     }
1446
1447     return content;
1448   }
1449
1450   bool CField::hasExpression(void) const
1451   {
1452     return (!expr.isEmpty() || !content.empty());
1453   }
1454
1455   DEFINE_REF_FUNC(Field,field)
1456} // namespace xios
Note: See TracBrowser for help on using the repository browser.