Index: /XIOS/dev/dev_olga/src/distribution_client.cpp
===================================================================
--- /XIOS/dev/dev_olga/src/distribution_client.cpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/distribution_client.cpp	(revision 1568)
@@ -36,5 +36,5 @@
   GlobalLocalMap void2 ;
   std::vector<int> void3 ;
-  std::vector<int> void4 ;
+  std::vector<bool> void4 ;
 
   globalLocalDataSendToServerMap_.swap(void1) ;
@@ -61,5 +61,5 @@
   // Then check mask of grid
   int gridDim = domList.size() * 2 + axisList.size();
-  grid->checkMask();
+//  grid->checkMask();
   switch (gridDim) {
     case 0:
@@ -282,6 +282,6 @@
 
         if ((iIdx >= nBeginLocal_[indexMap_[i]]) && (iIdx < nLocal_[indexMap_[i]]) &&
-           (jIdx >= nBeginLocal_[indexMap_[i]+1]) && (jIdx < nLocal_[indexMap_[i]+1]) &&
-           (domainMasks_[idxDomain](iIdx + jIdx*nLocal_[indexMap_[i]])))
+           (jIdx >= nBeginLocal_[indexMap_[i]+1]) && (jIdx < nLocal_[indexMap_[i]+1])
+           && (domainMasks_[idxDomain](iIdx + jIdx*nLocal_[indexMap_[i]])))
         {
           ++count;
@@ -325,5 +325,5 @@
       elementIndexData_[i].resize(dataNIndex_[i]);
       elementIndexData_[i] = false;
-      int iIdx = 0, count = 0, localIndex = 0;
+      int iIdx = 0, count = 0;
       for (int j = 0; j < dataNIndex_[i]; ++j)
       {
@@ -438,34 +438,20 @@
       currentIndex[0] = elementLocalIndex_[0](i);
 
-      // If defined, iterate on grid mask
-      if (!gridMask_.isEmpty())
-      {
-        for (int k = 0; k < this->numElement_; ++k)
-        {
-          gridMaskIndex += (currentIndex[k])*elementNLocal_[k];
-        }
-        if (gridMask_(gridMaskIndex)) ++indexLocalDataOnClientCount;
-      }
-      // If grid mask is not defined, iterate on elements' mask
-      else
-      {
-        bool maskTmp = true;
-        int idxDomain = 0, idxAxis = 0;
-        for (int elem = 0; elem < numElement_; ++elem)
-        {
-          if (2 == axisDomainOrder_(elem))
-          {
-            maskTmp = maskTmp && domainMasks_[idxDomain](currentIndex[elem]);
-            ++idxDomain;
-          }
-          else if (1 == axisDomainOrder_(elem))
-          {
-            maskTmp = maskTmp && axisMasks_[idxAxis](currentIndex[elem]);
-            ++idxAxis;
-          }
-        }
-        if (maskTmp) ++indexLocalDataOnClientCount;
-      }
-
+      bool maskTmp = true;
+      int idxDomain = 0, idxAxis = 0;
+      for (int elem = 0; elem < numElement_; ++elem)
+      {
+        if (2 == axisDomainOrder_(elem))
+        {
+          maskTmp = maskTmp && domainMasks_[idxDomain](currentIndex[elem]);
+          ++idxDomain;
+        }
+        else if (1 == axisDomainOrder_(elem))
+        {
+          maskTmp = maskTmp && axisMasks_[idxAxis](currentIndex[elem]);
+          ++idxAxis;
+        }
+      }
+      if (maskTmp) ++indexLocalDataOnClientCount;
     }
     idxLoop[0] += innerLoopSize;
@@ -535,33 +521,32 @@
         if (isCurrentIndexDataCorrect)
         {
-          int gridMaskIndex = 0;
-          for (int k = 0; k < this->numElement_; ++k)
+          bool maskTmp = true;
+          bool maskGridTmp = true;
+
+          // Domain or axis mask: only unmasked values will be sent
+          int idxDomain = 0, idxAxis = 0;
+          for (int elem = 0; elem < numElement_; ++elem)
           {
-            gridMaskIndex += (currentIndex[k])*elementNLocal_[k];
+            if (2 == axisDomainOrder_(elem))
+            {
+              maskTmp = maskTmp && domainMasks_[idxDomain](currentIndex[elem]);
+              ++idxDomain;
+            }
+            else if (1 == axisDomainOrder_(elem))
+            {
+              maskTmp = maskTmp && axisMasks_[idxAxis](currentIndex[elem]);
+              ++idxAxis;
+            }
           }
 
-          bool maskTmp = true;
-          // If defined, apply grid mask
-         if (!gridMask_.isEmpty())
+          // Grid mask: unmasked values will be replaces by NaN and then all values will be sent
+          if (!gridMask_.isEmpty())
           {
-            maskTmp =  gridMask_(gridMaskIndex);
-          }
-          // If grid mask is not defined, apply elements' mask
-          else
-          {
-            int idxDomain = 0, idxAxis = 0;
-            for (int elem = 0; elem < numElement_; ++elem)
+            int gridMaskIndex = 0;
+            for (int k = 0; k < this->numElement_; ++k)
             {
-              if (2 == axisDomainOrder_(elem))
-              {
-                maskTmp = maskTmp && domainMasks_[idxDomain](currentIndex[elem]);
-                ++idxDomain;
-              }
-              else if (1 == axisDomainOrder_(elem))
-              {
-                maskTmp = maskTmp && axisMasks_[idxAxis](currentIndex[elem]);
-                ++idxAxis;
-              }
+              gridMaskIndex += (currentIndex[k])*elementNLocal_[k];
             }
+            maskGridTmp =  gridMask_(gridMaskIndex);
           }
 
@@ -575,6 +560,9 @@
             globalDataIndex_[globalIndex] = indexLocalDataOnClientCount;
             localDataIndex_[indexLocalDataOnClientCount] = countLocalData;
+            if (maskGridTmp)
+              localMaskIndex_[indexLocalDataOnClientCount] = true;
+            else
+              localMaskIndex_[indexLocalDataOnClientCount] = false;
             globalLocalDataSendToServerMap_[globalIndex] = indexLocalDataOnClientCount;
-            localMaskIndex_[indexLocalDataOnClientCount] = gridMaskIndex;
             localMaskedDataIndex_[indexLocalDataOnClientCount] = indexLocalDataOnClientCount;
             ++indexLocalDataOnClientCount;
@@ -641,5 +629,5 @@
     return -1;
   }
-  int tempI = dataIndex + dataBegin;
+  int tempI = dataIndex; //+ dataBegin;
   if ((tempI < 0) || (tempI > ni))
     return -1;
@@ -675,5 +663,5 @@
   Return local mask index of client
 */
-const std::vector<int>& CDistributionClient::getLocalMaskIndexOnClient()
+const std::vector<bool>& CDistributionClient::getLocalMaskIndexOnClient()
 {
   if (!isComputed_) createGlobalIndexSendToServer();
@@ -691,2 +679,3 @@
 
 } // namespace xios
+
Index: /XIOS/dev/dev_olga/src/distribution_client.hpp
===================================================================
--- /XIOS/dev/dev_olga/src/distribution_client.hpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/distribution_client.hpp	(revision 1568)
@@ -33,5 +33,5 @@
 
   public:
-    /** Default constructor */        
+    /** Default constructor */
     CDistributionClient(int rank, CGrid* grid);
 
@@ -44,5 +44,5 @@
     GlobalLocalDataMap& getGlobalLocalDataSendToServer();
     GlobalLocalDataMap& getGlobalDataIndexOnClient();
-    const std::vector<int>& getLocalMaskIndexOnClient();
+    const std::vector<bool>& getLocalMaskIndexOnClient();
     const std::vector<int>& getLocalMaskedDataIndexOnClient();
 
@@ -83,6 +83,19 @@
     GlobalLocalDataMap globalLocalDataSendToServerMap_;
     GlobalLocalDataMap globalDataIndex_;
+
+    /*! Array holding masked data indexes.
+     * It includes:
+     *  masking on data (data_i/j_index or data_ni/nj and data_ibegin)
+     *  masking on grid elements (domains, axes or scalars)
+     * It DOES NOT include grid mask.
+     * The array size defines the data size entering the workflow. It is used by source filter of client or server1.
+    */
     std::vector<int> localDataIndex_;
-    std::vector<int> localMaskIndex_;
+
+    /*! Array holding grid mask.
+     * It is used by source filter of client for replacing unmasked data by NaN.
+    */
+    std::vector<bool> localMaskIndex_;
+
     std::vector<int> localMaskedDataIndex_;
 
@@ -172,2 +185,3 @@
 } // namespace xios
 #endif // __XIOS_DISTRIBUTIONCLIENT_HPP__
+
Index: /XIOS/dev/dev_olga/src/filter/source_filter.cpp
===================================================================
--- /XIOS/dev/dev_olga/src/filter/source_filter.cpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/filter/source_filter.cpp	(revision 1568)
@@ -7,5 +7,6 @@
 namespace xios
 {
-  CSourceFilter::CSourceFilter(CGarbageCollector& gc, CGrid* grid, bool compression, 
+  CSourceFilter::CSourceFilter(CGarbageCollector& gc, CGrid* grid,
+                               bool compression /*= true*/, bool mask /*= false*/,
                                const CDuration offset /*= NoneDu*/, bool manualTrigger /*= false*/,
                                bool hasMissingValue /*= false*/,
@@ -14,4 +15,5 @@
     , grid(grid)
     , compression(compression)
+    , mask(mask)
     , offset(offset)
     , hasMissingValue(hasMissingValue), defaultValue(defaultValue)
@@ -40,16 +42,10 @@
     }
     else
-      grid->inputField(data, packet->data);
-
-    
-    
-    // if (compression) grid->inputField(data, packet->data) ;
-    // else
-    // {
-    //   // just make a flat copy
-    //   CArray<double, N> data_tmp(data.copy()) ; // supress const attribute
-    //   CArray<double,1> dataTmp2(data_tmp.dataFirst(),shape(data.numElements()),neverDeleteData) ;
-    //   packet->data = dataTmp2 ;
-    // }
+    {
+      if (mask)
+        grid->maskField(data, packet->data);
+      else
+        grid->inputField(data, packet->data);
+    }
     // Convert missing values to NaN
     if (hasMissingValue)
Index: /XIOS/dev/dev_olga/src/filter/source_filter.hpp
===================================================================
--- /XIOS/dev/dev_olga/src/filter/source_filter.hpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/filter/source_filter.hpp	(revision 1568)
@@ -21,4 +21,6 @@
        * \param gc the garbage collector associated with this filter
        * \param grid the grid to which the data is attached
+       * \param compression
+       * \param mask
        * \param offset the offset applied to the timestamp of all packets
        * \param manualTrigger whether the output should be triggered manually
@@ -27,5 +29,6 @@
        */
       CSourceFilter(CGarbageCollector& gc, CGrid* grid,
-                    bool compression=true,
+                    bool compression = true,
+                    bool mask = false,
                     const CDuration offset = NoneDu, bool manualTrigger = false,
                     bool hasMissingValue = false,
@@ -61,9 +64,10 @@
 
     private:
-      CGrid* grid; //!< The grid attached to the data the filter can accept
-      const CDuration offset; //!< The offset applied to the timestamp of all packets
+      CGrid* grid;             //!< The grid attached to the data the filter can accept
+      const CDuration offset;  //!< The offset applied to the timestamp of all packets
       const bool hasMissingValue;
       const double defaultValue;
-      const bool compression ; //!< indicate if the data need to be compressed : on client size : true, on server side : false
+      const bool compression ; //!< indicates if data need to be compressed : on client side : true, on server side : false
+      const bool mask ;        //!< indicates whether grid mask should be applied (true for clients, false for servers)
   }; // class CSourceFilter
 } // namespace xios
Index: /XIOS/dev/dev_olga/src/node/axis.cpp
===================================================================
--- /XIOS/dev/dev_olga/src/node/axis.cpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/node/axis.cpp	(revision 1568)
@@ -285,12 +285,11 @@
       }
 
-      // Remove this check because it doen't make sense in case of a hole or overlapping axes
       if (!this->value.isEmpty())
       {
-//        StdSize true_size = value.numElements();
-//        if (this->n.getValue() != true_size)
-//          ERROR("CAxis::checkAttributes(void)",
-//                << "[ id = '" << getId() << "' , context = '" << CObjectFactory::GetCurrentContextId() << "' ] "
-//                << "The axis is wrongly defined, attribute 'value' has a different size (" << true_size << ") than the one defined by the \'size\' attribute (" << n.getValue() << ").");
+        StdSize true_size = value.numElements();
+        if (this->n.getValue() != true_size)
+          ERROR("CAxis::checkAttributes(void)",
+                << "[ id = '" << getId() << "' , context = '" << CObjectFactory::GetCurrentContextId() << "' ] "
+                << "The axis is wrongly defined, attribute 'value' has a different size (" << true_size << ") than the one defined by the \'size\' attribute (" << n.getValue() << ").");
         this->hasValue = true;
       }
@@ -310,4 +309,5 @@
       Check the validity of data and fill in values if any.
    */
+
    void CAxis::checkData()
    {
@@ -328,6 +328,30 @@
       {
         data_index.resize(data_n);
-        for (int i = 0; i < data_n; ++i) data_index(i) = i;
-      }
+        for (int i = 0; i < data_n; ++i)
+        {
+          if ((i+data_begin) >= 0 && (i+data_begin<n))
+            data_index(i) = i+data_begin;
+          else
+            data_index(i) = -1;
+        }
+      }
+      else
+      {
+        if (data_index.numElements() != data_n)
+        {
+          ERROR("CAxis::checkData(void)",
+                << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
+                << "The size of data_index = "<< data_index.numElements() << "is not equal to the data size data_n = " << data_n.getValue() << ").");
+        }
+        data_index.resize(data_n);
+        for (int i = 0; i < data_n; ++i)
+        {
+          if ((i+data_begin) >= 0 && (i+data_begin<n))
+            data_index(i) = data_index(i);
+          else
+            data_index(i) = -1;
+        }
+      }
+
    }
 
@@ -344,18 +368,17 @@
       if (!mask.isEmpty())
       {
-         if (mask.extent(0) != n)
-           ERROR("CAxis::checkMask(void)",
-                 << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
-                 << "The mask does not have the same size as the local domain." << std::endl
-                 << "Local size is " << n.getValue() << "." << std::endl
-                 << "Mask size is " << mask.extent(0) << ".");
-      }
-      else // (mask.isEmpty())
-      { // If no mask was defined, we create a default one without any masked point.
-         mask.resize(n);
-         for (int i = 0; i < n; ++i)
-         {
-           mask(i) = true;
-         }
+        if (mask.extent(0) != n)
+        {
+          ERROR("CAxis::checkMask(void)",
+              << "[ id = " << this->getId() << " , context = '" << CObjectFactory::GetCurrentContextId() << " ] "
+              << "The mask does not have the same size as the local domain." << std::endl
+              << "Local size is " << n.getValue() << "." << std::endl
+              << "Mask size is " << mask.extent(0) << ".");
+        }
+      }
+      else
+      {
+        mask.resize(n);
+        mask = true;
       }
    }
@@ -547,18 +570,18 @@
 
         // Calculate the compressed index if any
-        std::set<int> writtenInd;
-        if (isCompressible_)
-        {
-          for (int idx = 0; idx < data_index.numElements(); ++idx)
-          {
-            int ind = CDistributionClient::getAxisIndex(data_index(idx), data_begin, ni);
-
-            if (ind >= 0 && ind < ni && mask(ind))
-            {
-              ind += ibegin;
-              writtenInd.insert(ind);
-            }
-          }
-        }
+//        std::set<int> writtenInd;
+//        if (isCompressible_)
+//        {
+//          for (int idx = 0; idx < data_index.numElements(); ++idx)
+//          {
+//            int ind = CDistributionClient::getAxisIndex(data_index(idx), data_begin, ni);
+//
+//            if (ind >= 0 && ind < ni && mask(ind))
+//            {
+//              ind += ibegin;
+//              writtenInd.insert(ind);
+//            }
+//          }
+//        }
 
         // Compute the global index of the current client (process) hold
@@ -725,4 +748,24 @@
         }                 
       }
+//
+//      nbWritten = 0;
+//      for (int idx = 0; idx < data_index.numElements(); ++idx)
+//      {
+//        if (localGlobalIndexMap.end() != localGlobalIndexMap.find(data_index(idx)))
+//        {
+//          ++nbWritten;
+//        }
+//      }
+//
+//      compressedIndexToWriteOnServer[writtenCommSize].resize(nbWritten);
+//      nbWritten = 0;
+//      for (int idx = 0; idx < data_index.numElements(); ++idx)
+//      {
+//        if (localGlobalIndexMap.end() != localGlobalIndexMap.find(data_index(idx)))
+//        {
+//          compressedIndexToWriteOnServer[writtenCommSize](nbWritten) = localGlobalIndexMap[data_index(idx)];
+//          ++nbWritten;
+//        }
+//      }
 
       nbWritten = 0;
@@ -959,11 +1002,11 @@
 
   /*
-    Send attributes of axis from a group of client to other group of clients/servers 
-    on supposing that these attributes are distributed among the clients of the sending group
-    In the future, if new attributes are added, they should also be processed in this function
+    Send axis attributes from a group of clients to another group of clients/servers
+    supposing that these attributes are distributed among the clients of the sending group
+    In future, if new attributes are added, they should also be processed in this function
   */
   void CAxis::sendDistributedAttributes(void)
   {
-    int ns, n, i, j, ind, nv, idx;
+    int ind, idx;
     std::list<CContextClient*>::iterator it;
 
@@ -982,4 +1025,5 @@
       list<CArray<string,1> > list_label;
 
+      // Cut off the ghost points
       int nbIndex = index.numElements();
       CArray<int,1> dataIndex(nbIndex);
@@ -988,5 +1032,5 @@
       {
         if (0 <= data_index(idx) && data_index(idx) < nbIndex)
-          dataIndex(idx) = 1;
+          dataIndex(data_index(idx)) = 1;
       }
 
@@ -995,5 +1039,5 @@
       for (int k = 0; k < connectedServerRank_[nbServer].size(); ++k)
       {
-        int nbData = 0;
+        int nbData = 0, nbDataCount = 0;
         int rank = connectedServerRank_[nbServer][k];
         it = indSrv_[nbServer].find(rank);
@@ -1002,5 +1046,5 @@
 
         list_indi.push_back(CArray<int,1>(nbData));
-        list_dataInd.push_back(CArray<int,1>(nbData));        
+        list_dataInd.push_back(CArray<int,1>(nbData));
         list_mask.push_back(CArray<bool,1>(nbData));
 
@@ -1015,8 +1059,9 @@
 
         CArray<int,1>& indi = list_indi.back();
-        CArray<int,1>& dataIndi = list_dataInd.back();        
+        CArray<int,1>& dataIndi = list_dataInd.back();
+        dataIndi = -1;
         CArray<bool,1>& maskIndi = list_mask.back();
 
-        for (n = 0; n < nbData; ++n)
+        for (int n = 0; n < nbData; ++n)
         {
           idx = static_cast<int>(it->second[n]);
@@ -1162,4 +1207,5 @@
     nonCompressedData = -1;   
     mask.resize(nbData);
+    mask = true;
     if (hasValue)
       value.resize(nbData);
@@ -1199,9 +1245,9 @@
     }
     
-    int nbCompressedData = 0; 
+    int nbCompressedData = 0;
     for (idx = 0; idx < nonCompressedData.numElements(); ++idx)
     {
       if (0 <= nonCompressedData(idx))
-        ++nbCompressedData;        
+        ++nbCompressedData;
     }
 
@@ -1213,9 +1259,10 @@
       {
         data_index(nbCompressedData) = idx % n;
-        ++nbCompressedData;        
+        ++nbCompressedData;
       }
     }
 
     data_begin.setValue(0);
+    data_n.setValue(data_index.numElements());
   }
 
@@ -1330,5 +1377,5 @@
       clientsSet.insert(contextClient);
     }
-}
+  }
 
   void CAxis::parse(xml::CXMLNode & node)
Index: /XIOS/dev/dev_olga/src/node/axis.hpp
===================================================================
--- /XIOS/dev/dev_olga/src/node/axis.hpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/node/axis.hpp	(revision 1568)
@@ -126,5 +126,5 @@
         bool hasLabel;
 
-        CArray<size_t,1> localIndexToWriteOnServer;        
+        CArray<size_t,1> localIndexToWriteOnServer;
 
       private:
Index: /XIOS/dev/dev_olga/src/node/domain.cpp
===================================================================
--- /XIOS/dev/dev_olga/src/node/domain.cpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/node/domain.cpp	(revision 1568)
@@ -2447,7 +2447,7 @@
          {
            iIndex = (index%ni_glo)-ibegin;
-           iIndex = (iIndex < 0) ? 0 : iIndex;
+           iIndex = (iIndex < 0) ? 0 : iIndex;    // ??
            jIndex = (index/ni_glo)-jbegin;
-           jIndex = (jIndex < 0) ? 0 : jIndex;
+           jIndex = (jIndex < 0) ? 0 : jIndex;    // ??
            nbIndLoc = iIndex + ni * jIndex;
            if (nbIndLoc < nbIndexGlobMax)
@@ -2896,9 +2896,9 @@
       buffer >> recvDataJIndex[i];
     }
-   
+
     int nbIndex = i_index.numElements();
     CArray<int,1> dataIIndex(nbIndex), dataJIndex(nbIndex);
     dataIIndex = -1; dataJIndex = -1;
-     
+
     nbIndex = 0;
     for (i = 0; i < nbReceived; ++i)
Index: /XIOS/dev/dev_olga/src/node/field.cpp
===================================================================
--- /XIOS/dev/dev_olga/src/node/field.cpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/node/field.cpp	(revision 1568)
@@ -1027,5 +1027,5 @@
      {
         if (!instantDataFilter)
-          instantDataFilter = clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid,true));
+          instantDataFilter = clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, false));
 
 
@@ -1043,5 +1043,5 @@
      {
        if (!instantDataFilter)
-         instantDataFilter = clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true));
+         instantDataFilter = clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, false, false));
 
              // If the field data is to be read by the client or/and written to a file
@@ -1089,5 +1089,5 @@
          {
            checkTimeAttributes();
-           instantDataFilter = serverSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, freq_offset, true,
+           instantDataFilter = serverSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, false, freq_offset, true,
                                                                                                        detectMissingValues, defaultValue));
          }
@@ -1095,5 +1095,5 @@
          {
             if (check_if_active.isEmpty()) check_if_active = false; 
-            instantDataFilter = clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, false, NoneDu, false,
+            instantDataFilter = clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, false, true, NoneDu, false,
                                                                                                         detectMissingValues, defaultValue));
          }
@@ -1176,5 +1176,5 @@
          {
            checkTimeAttributes();
-           serverSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, freq_offset, true,
+           serverSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, false, freq_offset, true,
                                                                                    detectMissingValues, defaultValue));
          }
@@ -1193,5 +1193,5 @@
          {
            if (check_if_active.isEmpty()) check_if_active = false;
-           clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, NoneDu, false,
+           clientSourceFilter = std::shared_ptr<CSourceFilter>(new CSourceFilter(gc, grid, true, true, NoneDu, false,
                                                                                    detectMissingValues, defaultValue));
          }
Index: /XIOS/dev/dev_olga/src/node/grid.cpp
===================================================================
--- /XIOS/dev/dev_olga/src/node/grid.cpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/node/grid.cpp	(revision 1568)
@@ -698,7 +698,12 @@
      clientDistribution_ = new CDistributionClient(rank, this);
      // Get local data index on client
-     storeIndex_client.resize(clientDistribution_->getLocalDataIndexOnClient().size());
-     int nbStoreIndex = storeIndex_client.numElements();
-     for (int idx = 0; idx < nbStoreIndex; ++idx) storeIndex_client(idx) = (clientDistribution_->getLocalDataIndexOnClient())[idx];
+     int nbStoreIndex = clientDistribution_->getLocalDataIndexOnClient().size();
+     storeIndex_client.resize(nbStoreIndex);
+     storeMask_client.resize(nbStoreIndex);
+     for (int idx = 0; idx < nbStoreIndex; ++idx)
+     {
+       storeIndex_client(idx) = (clientDistribution_->getLocalDataIndexOnClient())[idx];
+       storeMask_client(idx) = (clientDistribution_->getLocalMaskIndexOnClient())[idx];
+     }
 
      if (0 == serverDistribution_) isDataDistributed_= clientDistribution_->isDataDistributed();
@@ -1293,4 +1298,13 @@
    }
 
+   void CGrid::maskField_arr(const double* const data, CArray<double, 1>& stored) const
+   {
+      const StdSize size = storeIndex_client.numElements();
+      stored.resize(size);
+
+      const double nanValue = std::numeric_limits<double>::quiet_NaN();
+      for(StdSize i = 0; i < size; i++) stored(i) = (storeMask_client(i)) ? data[storeIndex_client(i)] : nanValue;
+   }
+
    void CGrid::uncompressField_arr(const double* const data, CArray<double, 1>& out) const
    {
@@ -1749,39 +1763,40 @@
         }
 
-        modifyMaskSize(nSize, false);
-
-        // These below codes are reserved for future
-        CDistributionServer srvDist(server->intraCommRank, nBegin, nSize, nBeginGlobal, nGlob); 
-        map<int, CArray<size_t, 1> >::iterator itb = outGlobalIndexFromClient.begin(),
-                                               ite = outGlobalIndexFromClient.end(), it;  
-        const CDistributionServer::GlobalLocalMap&  globalLocalMask = srvDist.getGlobalLocalIndex();
-        CDistributionServer::GlobalLocalMap::const_iterator itSrv;
-        size_t nb = 0;
-        for (it = itb; it != ite; ++it)
-        {
-          CArray<size_t,1>& globalInd = it->second;
-          for (size_t idx = 0; idx < globalInd.numElements(); ++idx)
-          {
-            if (globalLocalMask.end() != globalLocalMask.find(globalInd(idx))) ++nb;
-          }
-        }
-        
-        CArray<int,1> indexToModify(nb);
-        nb = 0;    
-        for (it = itb; it != ite; ++it)
-        {
-          CArray<size_t,1>& globalInd = it->second;
-          for (size_t idx = 0; idx < globalInd.numElements(); ++idx)
-          {
-            itSrv = globalLocalMask.find(globalInd(idx));
-            if (globalLocalMask.end() != itSrv) 
-            {
-              indexToModify(nb) = itSrv->second;
-              ++nb;
-            }
-          }
-        }
-
-        modifyMask(indexToModify, true);
+//        modifyMaskSize(nSize, false);
+//        modifyMaskSize(nSize, true);   // Grid mask should be gone beyond client source filter
+//
+//        // These below codes are reserved for future
+//        CDistributionServer srvDist(server->intraCommRank, nBegin, nSize, nBeginGlobal, nGlob);
+//        map<int, CArray<size_t, 1> >::iterator itb = outGlobalIndexFromClient.begin(),
+//                                               ite = outGlobalIndexFromClient.end(), it;
+//        const CDistributionServer::GlobalLocalMap&  globalLocalMask = srvDist.getGlobalLocalIndex();
+//        CDistributionServer::GlobalLocalMap::const_iterator itSrv;
+//        size_t nb = 0;
+//        for (it = itb; it != ite; ++it)
+//        {
+//          CArray<size_t,1>& globalInd = it->second;
+//          for (size_t idx = 0; idx < globalInd.numElements(); ++idx)
+//          {
+//            if (globalLocalMask.end() != globalLocalMask.find(globalInd(idx))) ++nb;
+//          }
+//        }
+//
+//        CArray<int,1> indexToModify(nb);
+//        nb = 0;
+//        for (it = itb; it != ite; ++it)
+//        {
+//          CArray<size_t,1>& globalInd = it->second;
+//          for (size_t idx = 0; idx < globalInd.numElements(); ++idx)
+//          {
+//            itSrv = globalLocalMask.find(globalInd(idx));
+//            if (globalLocalMask.end() != itSrv)
+//            {
+//              indexToModify(nb) = itSrv->second;
+//              ++nb;
+//            }
+//          }
+//        }
+
+//        modifyMask(indexToModify, true);
       }
 
Index: /XIOS/dev/dev_olga/src/node/grid.hpp
===================================================================
--- /XIOS/dev/dev_olga/src/node/grid.hpp	(revision 1567)
+++ /XIOS/dev/dev_olga/src/node/grid.hpp	(revision 1568)
@@ -94,4 +94,6 @@
          template <int n>
          void inputField(const CArray<double,n>& field, CArray<double,1>& stored) const;
+         template <int n>
+         void maskField(const CArray<double,n>& field, CArray<double,1>& stored) const;
          template <int n>
          void outputField(const CArray<double,1>& stored, CArray<double,n>& field) const;  
@@ -208,4 +210,5 @@
       public:
          CArray<int, 1> storeIndex_client;
+         CArray<bool, 1> storeMask_client;
 
 /** Map containing indexes that will be sent in sendIndex(). */
@@ -273,4 +276,5 @@
         void restoreField_arr(const CArray<double, 1>& stored, double* const data) const;
         void uncompressField_arr(const double* const data, CArray<double, 1>& outData) const;
+        void maskField_arr(const double* const data, CArray<double, 1>& stored) const;
 
         void setVirtualDomainGroup(CDomainGroup* newVDomainGroup);
@@ -385,4 +389,18 @@
 
    template <int n>
+   void CGrid::maskField(const CArray<double,n>& field, CArray<double,1>& stored) const
+   {
+//#ifdef __XIOS_DEBUG
+      if (this->getDataSize() != field.numElements())
+         ERROR("void CGrid::inputField(const  CArray<double,n>& field, CArray<double,1>& stored) const",
+                << "[ Awaiting data of size = " << this->getDataSize() << ", "
+                << "Received data size = "      << field.numElements() << " ] "
+                << "The data array does not have the right size! "
+                << "Grid = " << this->getId())
+//#endif
+      this->maskField_arr(field.dataFirst(), stored);
+   }
+
+   template <int n>
    void CGrid::outputField(const CArray<double,1>& stored, CArray<double,n>& field) const
    {
@@ -417,101 +435,97 @@
                              bool createMask)
    {
-     if (!gridMask.isEmpty() || createMask)
-     {
-       int idx = 0;
-       int numElement = axisDomainOrder.numElements();
-       int dim = domainMasks.size() * 2 + axisMasks.size();
-       std::vector<CDomain*> domainP = this->getDomains();
-       std::vector<CAxis*> axisP = this->getAxis();
-
-       std::vector<int> idxLoop(dim,0), indexMap(numElement), eachDimSize(dim);
-       std::vector<int> currentIndex(dim);
-       int idxDomain = 0, idxAxis = 0;
+     int idx = 0;
+     int numElement = axisDomainOrder.numElements();
+     int dim = domainMasks.size() * 2 + axisMasks.size();
+     std::vector<CDomain*> domainP = this->getDomains();
+     std::vector<CAxis*> axisP = this->getAxis();
+
+     std::vector<int> idxLoop(dim,0), indexMap(numElement), eachDimSize(dim);
+     std::vector<int> currentIndex(dim);
+     int idxDomain = 0, idxAxis = 0;
+    for (int i = 0; i < numElement; ++i)
+    {
+      indexMap[i] = idx;
+      if (2 == axisDomainOrder(i)) {
+          eachDimSize[indexMap[i]]   = domainP[idxDomain]->ni;
+          eachDimSize[indexMap[i]+1] = domainP[idxDomain]->nj;
+          idx += 2; ++idxDomain;
+      }
+      else if (1 == axisDomainOrder(i)) {
+//        eachDimSize[indexMap[i]] = axisMasks[idxAxis]->numElements();
+        eachDimSize[indexMap[i]] = axisP[idxAxis]->n;
+        ++idx; ++idxAxis;
+      }
+      else {};
+    }
+
+    if (!gridMask.isEmpty() && !createMask)
+    {
+      for (int i = 0; i < dim; ++i)
+      {
+        if (gridMask.extent(i) != eachDimSize[i])
+          ERROR("CGrid::checkMask(void)",
+                << "The mask has one dimension whose size is different from the one of the local grid." << std::endl
+                << "Local size of dimension " << i << " is " << eachDimSize[i] << "." << std::endl
+                << "Mask size for dimension " << i << " is " << gridMask.extent(i) << "." << std::endl
+                << "Grid = " << this->getId())
+      }
+    }
+    else {
+        CArrayBoolTraits<CArray<bool,N> >::resizeArray(gridMask,eachDimSize);
+        gridMask = true;
+    }
+
+    int ssize = gridMask.numElements();
+    idx = 0;
+    while (idx < ssize)
+    {
+      for (int i = 0; i < dim-1; ++i)
+      {
+        if (idxLoop[i] == eachDimSize[i])
+        {
+          idxLoop[i] = 0;
+          ++idxLoop[i+1];
+        }
+      }
+
+      // Find out outer index
+      idxDomain = idxAxis = 0;
+      bool maskValue = true;
       for (int i = 0; i < numElement; ++i)
       {
-        indexMap[i] = idx;
-        if (2 == axisDomainOrder(i)) {
-            eachDimSize[indexMap[i]]   = domainP[idxDomain]->ni;
-            eachDimSize[indexMap[i]+1] = domainP[idxDomain]->nj;
-            idx += 2; ++idxDomain;
+        if (2 == axisDomainOrder(i))
+        {
+          int idxTmp = idxLoop[indexMap[i]] + idxLoop[indexMap[i]+1] * eachDimSize[indexMap[i]];
+          if (idxTmp < (*domainMasks[idxDomain]).numElements())
+            maskValue = maskValue && (*domainMasks[idxDomain])(idxTmp);
+          else
+            maskValue = false;
+          ++idxDomain;
         }
-        else if (1 == axisDomainOrder(i)) {
-  //        eachDimSize[indexMap[i]] = axisMasks[idxAxis]->numElements();
-          eachDimSize[indexMap[i]] = axisP[idxAxis]->n;
-          ++idx; ++idxAxis;
+        else if (1 == axisDomainOrder(i))
+        {
+          int idxTmp = idxLoop[indexMap[i]];
+          if (idxTmp < (*axisMasks[idxAxis]).numElements())
+            maskValue = maskValue && (*axisMasks[idxAxis])(idxTmp);
+          else
+            maskValue = false;
+
+          ++idxAxis;
         }
-        else {};
-      }
-
-//      if (!gridMask.isEmpty() && !createMask)
-      if (!createMask)
+      }
+
+      int maskIndex = idxLoop[0];
+      int mulDim = 1;
+      for (int k = 1; k < dim; ++k)
       {
-        for (int i = 0; i < dim; ++i)
-        {
-          if (gridMask.extent(i) != eachDimSize[i])
-            ERROR("CGrid::checkMask(void)",
-                  << "The mask has one dimension whose size is different from the one of the local grid." << std::endl
-                  << "Local size of dimension " << i << " is " << eachDimSize[i] << "." << std::endl
-                  << "Mask size for dimension " << i << " is " << gridMask.extent(i) << "." << std::endl
-                  << "Grid = " << this->getId())
-        }
-      }
-      else {
-          CArrayBoolTraits<CArray<bool,N> >::resizeArray(gridMask,eachDimSize);
-          gridMask = true;
-      }
-
-      int ssize = gridMask.numElements();
-      idx = 0;
-      while (idx < ssize)
-      {
-        for (int i = 0; i < dim-1; ++i)
-        {
-          if (idxLoop[i] == eachDimSize[i])
-          {
-            idxLoop[i] = 0;
-            ++idxLoop[i+1];
-          }
-        }
-
-        // Find out outer index
-        idxDomain = idxAxis = 0;
-        bool maskValue = true;
-        for (int i = 0; i < numElement; ++i)
-        {
-          if (2 == axisDomainOrder(i))
-          {
-            int idxTmp = idxLoop[indexMap[i]] + idxLoop[indexMap[i]+1] * eachDimSize[indexMap[i]];
-            if (idxTmp < (*domainMasks[idxDomain]).numElements())
-              maskValue = maskValue && (*domainMasks[idxDomain])(idxTmp);
-            else
-              maskValue = false;
-            ++idxDomain;
-          }
-          else if (1 == axisDomainOrder(i))
-          {
-            int idxTmp = idxLoop[indexMap[i]];
-            if (idxTmp < (*axisMasks[idxAxis]).numElements())
-              maskValue = maskValue && (*axisMasks[idxAxis])(idxTmp);
-            else
-              maskValue = false;
-
-            ++idxAxis;
-          }
-        }
-
-        int maskIndex = idxLoop[0];
-        int mulDim = 1;
-        for (int k = 1; k < dim; ++k)
-        {
-          mulDim *= eachDimSize[k-1];
-          maskIndex += idxLoop[k]*mulDim;
-        }
-        *(gridMask.dataFirst()+maskIndex) &= maskValue;
-
-        ++idxLoop[0];
-        ++idx;
-      }
-     }
+        mulDim *= eachDimSize[k-1];
+        maskIndex += idxLoop[k]*mulDim;
+      }
+      *(gridMask.dataFirst()+maskIndex) &= maskValue;
+
+      ++idxLoop[0];
+      ++idx;
+    }
    }
 
