Ignore:
Timestamp:
08/28/18 16:38:30 (6 years ago)
Author:
oabramkina
Message:

Grid mask is now applied in the source filter of clients: unmasked values are replaced by NaN. It is not reconstructed any more by servers.

This needs to be tested more rigorously before commiting to trunk.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • XIOS/dev/dev_olga/src/node/grid.hpp

    r1564 r1568  
    9494         template <int n> 
    9595         void inputField(const CArray<double,n>& field, CArray<double,1>& stored) const; 
     96         template <int n> 
     97         void maskField(const CArray<double,n>& field, CArray<double,1>& stored) const; 
    9698         template <int n> 
    9799         void outputField(const CArray<double,1>& stored, CArray<double,n>& field) const;   
     
    208210      public: 
    209211         CArray<int, 1> storeIndex_client; 
     212         CArray<bool, 1> storeMask_client; 
    210213 
    211214/** Map containing indexes that will be sent in sendIndex(). */ 
     
    273276        void restoreField_arr(const CArray<double, 1>& stored, double* const data) const; 
    274277        void uncompressField_arr(const double* const data, CArray<double, 1>& outData) const; 
     278        void maskField_arr(const double* const data, CArray<double, 1>& stored) const; 
    275279 
    276280        void setVirtualDomainGroup(CDomainGroup* newVDomainGroup); 
     
    385389 
    386390   template <int n> 
     391   void CGrid::maskField(const CArray<double,n>& field, CArray<double,1>& stored) const 
     392   { 
     393//#ifdef __XIOS_DEBUG 
     394      if (this->getDataSize() != field.numElements()) 
     395         ERROR("void CGrid::inputField(const  CArray<double,n>& field, CArray<double,1>& stored) const", 
     396                << "[ Awaiting data of size = " << this->getDataSize() << ", " 
     397                << "Received data size = "      << field.numElements() << " ] " 
     398                << "The data array does not have the right size! " 
     399                << "Grid = " << this->getId()) 
     400//#endif 
     401      this->maskField_arr(field.dataFirst(), stored); 
     402   } 
     403 
     404   template <int n> 
    387405   void CGrid::outputField(const CArray<double,1>& stored, CArray<double,n>& field) const 
    388406   { 
     
    417435                             bool createMask) 
    418436   { 
    419      if (!gridMask.isEmpty() || createMask) 
    420      { 
    421        int idx = 0; 
    422        int numElement = axisDomainOrder.numElements(); 
    423        int dim = domainMasks.size() * 2 + axisMasks.size(); 
    424        std::vector<CDomain*> domainP = this->getDomains(); 
    425        std::vector<CAxis*> axisP = this->getAxis(); 
    426  
    427        std::vector<int> idxLoop(dim,0), indexMap(numElement), eachDimSize(dim); 
    428        std::vector<int> currentIndex(dim); 
    429        int idxDomain = 0, idxAxis = 0; 
     437     int idx = 0; 
     438     int numElement = axisDomainOrder.numElements(); 
     439     int dim = domainMasks.size() * 2 + axisMasks.size(); 
     440     std::vector<CDomain*> domainP = this->getDomains(); 
     441     std::vector<CAxis*> axisP = this->getAxis(); 
     442 
     443     std::vector<int> idxLoop(dim,0), indexMap(numElement), eachDimSize(dim); 
     444     std::vector<int> currentIndex(dim); 
     445     int idxDomain = 0, idxAxis = 0; 
     446    for (int i = 0; i < numElement; ++i) 
     447    { 
     448      indexMap[i] = idx; 
     449      if (2 == axisDomainOrder(i)) { 
     450          eachDimSize[indexMap[i]]   = domainP[idxDomain]->ni; 
     451          eachDimSize[indexMap[i]+1] = domainP[idxDomain]->nj; 
     452          idx += 2; ++idxDomain; 
     453      } 
     454      else if (1 == axisDomainOrder(i)) { 
     455//        eachDimSize[indexMap[i]] = axisMasks[idxAxis]->numElements(); 
     456        eachDimSize[indexMap[i]] = axisP[idxAxis]->n; 
     457        ++idx; ++idxAxis; 
     458      } 
     459      else {}; 
     460    } 
     461 
     462    if (!gridMask.isEmpty() && !createMask) 
     463    { 
     464      for (int i = 0; i < dim; ++i) 
     465      { 
     466        if (gridMask.extent(i) != eachDimSize[i]) 
     467          ERROR("CGrid::checkMask(void)", 
     468                << "The mask has one dimension whose size is different from the one of the local grid." << std::endl 
     469                << "Local size of dimension " << i << " is " << eachDimSize[i] << "." << std::endl 
     470                << "Mask size for dimension " << i << " is " << gridMask.extent(i) << "." << std::endl 
     471                << "Grid = " << this->getId()) 
     472      } 
     473    } 
     474    else { 
     475        CArrayBoolTraits<CArray<bool,N> >::resizeArray(gridMask,eachDimSize); 
     476        gridMask = true; 
     477    } 
     478 
     479    int ssize = gridMask.numElements(); 
     480    idx = 0; 
     481    while (idx < ssize) 
     482    { 
     483      for (int i = 0; i < dim-1; ++i) 
     484      { 
     485        if (idxLoop[i] == eachDimSize[i]) 
     486        { 
     487          idxLoop[i] = 0; 
     488          ++idxLoop[i+1]; 
     489        } 
     490      } 
     491 
     492      // Find out outer index 
     493      idxDomain = idxAxis = 0; 
     494      bool maskValue = true; 
    430495      for (int i = 0; i < numElement; ++i) 
    431496      { 
    432         indexMap[i] = idx; 
    433         if (2 == axisDomainOrder(i)) { 
    434             eachDimSize[indexMap[i]]   = domainP[idxDomain]->ni; 
    435             eachDimSize[indexMap[i]+1] = domainP[idxDomain]->nj; 
    436             idx += 2; ++idxDomain; 
     497        if (2 == axisDomainOrder(i)) 
     498        { 
     499          int idxTmp = idxLoop[indexMap[i]] + idxLoop[indexMap[i]+1] * eachDimSize[indexMap[i]]; 
     500          if (idxTmp < (*domainMasks[idxDomain]).numElements()) 
     501            maskValue = maskValue && (*domainMasks[idxDomain])(idxTmp); 
     502          else 
     503            maskValue = false; 
     504          ++idxDomain; 
    437505        } 
    438         else if (1 == axisDomainOrder(i)) { 
    439   //        eachDimSize[indexMap[i]] = axisMasks[idxAxis]->numElements(); 
    440           eachDimSize[indexMap[i]] = axisP[idxAxis]->n; 
    441           ++idx; ++idxAxis; 
     506        else if (1 == axisDomainOrder(i)) 
     507        { 
     508          int idxTmp = idxLoop[indexMap[i]]; 
     509          if (idxTmp < (*axisMasks[idxAxis]).numElements()) 
     510            maskValue = maskValue && (*axisMasks[idxAxis])(idxTmp); 
     511          else 
     512            maskValue = false; 
     513 
     514          ++idxAxis; 
    442515        } 
    443         else {}; 
    444       } 
    445  
    446 //      if (!gridMask.isEmpty() && !createMask) 
    447       if (!createMask) 
     516      } 
     517 
     518      int maskIndex = idxLoop[0]; 
     519      int mulDim = 1; 
     520      for (int k = 1; k < dim; ++k) 
    448521      { 
    449         for (int i = 0; i < dim; ++i) 
    450         { 
    451           if (gridMask.extent(i) != eachDimSize[i]) 
    452             ERROR("CGrid::checkMask(void)", 
    453                   << "The mask has one dimension whose size is different from the one of the local grid." << std::endl 
    454                   << "Local size of dimension " << i << " is " << eachDimSize[i] << "." << std::endl 
    455                   << "Mask size for dimension " << i << " is " << gridMask.extent(i) << "." << std::endl 
    456                   << "Grid = " << this->getId()) 
    457         } 
    458       } 
    459       else { 
    460           CArrayBoolTraits<CArray<bool,N> >::resizeArray(gridMask,eachDimSize); 
    461           gridMask = true; 
    462       } 
    463  
    464       int ssize = gridMask.numElements(); 
    465       idx = 0; 
    466       while (idx < ssize) 
    467       { 
    468         for (int i = 0; i < dim-1; ++i) 
    469         { 
    470           if (idxLoop[i] == eachDimSize[i]) 
    471           { 
    472             idxLoop[i] = 0; 
    473             ++idxLoop[i+1]; 
    474           } 
    475         } 
    476  
    477         // Find out outer index 
    478         idxDomain = idxAxis = 0; 
    479         bool maskValue = true; 
    480         for (int i = 0; i < numElement; ++i) 
    481         { 
    482           if (2 == axisDomainOrder(i)) 
    483           { 
    484             int idxTmp = idxLoop[indexMap[i]] + idxLoop[indexMap[i]+1] * eachDimSize[indexMap[i]]; 
    485             if (idxTmp < (*domainMasks[idxDomain]).numElements()) 
    486               maskValue = maskValue && (*domainMasks[idxDomain])(idxTmp); 
    487             else 
    488               maskValue = false; 
    489             ++idxDomain; 
    490           } 
    491           else if (1 == axisDomainOrder(i)) 
    492           { 
    493             int idxTmp = idxLoop[indexMap[i]]; 
    494             if (idxTmp < (*axisMasks[idxAxis]).numElements()) 
    495               maskValue = maskValue && (*axisMasks[idxAxis])(idxTmp); 
    496             else 
    497               maskValue = false; 
    498  
    499             ++idxAxis; 
    500           } 
    501         } 
    502  
    503         int maskIndex = idxLoop[0]; 
    504         int mulDim = 1; 
    505         for (int k = 1; k < dim; ++k) 
    506         { 
    507           mulDim *= eachDimSize[k-1]; 
    508           maskIndex += idxLoop[k]*mulDim; 
    509         } 
    510         *(gridMask.dataFirst()+maskIndex) &= maskValue; 
    511  
    512         ++idxLoop[0]; 
    513         ++idx; 
    514       } 
    515      } 
     522        mulDim *= eachDimSize[k-1]; 
     523        maskIndex += idxLoop[k]*mulDim; 
     524      } 
     525      *(gridMask.dataFirst()+maskIndex) &= maskValue; 
     526 
     527      ++idxLoop[0]; 
     528      ++idx; 
     529    } 
    516530   } 
    517531 
Note: See TracChangeset for help on using the changeset viewer.