Changeset 1480 for XIOS/dev

Timestamp:

04/13/18 10:33:54 (6 years ago)

Author:

ymipsl

Message:

Fix for transformation & interpolation : better detection of missing values.

YM

Location:

XIOS/dev/XIOS_DEV_CMIP6/src/transformation

Files:

• domain_algorithm_interpolate.cpp (modified) (3 diffs)
• domain_algorithm_interpolate.hpp (modified) (1 diff)
• generic_algorithm_transformation.cpp (modified) (1 diff)

XIOS/dev/XIOS_DEV_CMIP6/src/transformation/domain_algorithm_interpolate.cpp

@@ -912 +912 @@
      }
 
+    if (firstPass)
+    {
+      allMissing.resize(dataOut.numElements()) ;
+      allMissing=true ;
+    }
+
     for (int idx = 0; idx < nbLocalIndex; ++idx)
     {
       if (NumTraits<double>::isNan(*(dataInput + idx)))
       {
-        flagInitial[localIndex[idx].first] = false;
+        allMissing(localIndex[idx].first) = allMissing(localIndex[idx].first) && true;
         if (renormalize) renormalizationFactor(localIndex[idx].first)-=localIndex[idx].second ;
       }
@@ -922 +928 @@
       {
         dataOut(localIndex[idx].first) += *(dataInput + idx) * localIndex[idx].second;
-        flagInitial[localIndex[idx].first] = true; // Reset flag to indicate not all data source are nan
+        allMissing(localIndex[idx].first) = allMissing(localIndex[idx].first) && false; // Reset flag to indicate not all data source are nan
       }
     }
 
-    // If all data source are nan then data destination must be nan
+  }
+  else
+  {
     for (int idx = 0; idx < nbLocalIndex; ++idx)
     {
-      if (!flagInitial[localIndex[idx].first])
-        dataOut(localIndex[idx].first) = defaultValue;
-    }
-  }
-  else
-  {
-    for (int idx = 0; idx < nbLocalIndex; ++idx)
-    {
       dataOut(localIndex[idx].first) += *(dataInput + idx) * localIndex[idx].second;
     }
@@ -944 +944 @@
 void CDomainAlgorithmInterpolate::updateData(CArray<double,1>& dataOut)
 {
-  if (detectMissingValue && renormalize)
-  {
-    if (renormalizationFactor.numElements()>0) dataOut/=renormalizationFactor ; // In some case, process doesn't received any data for interpolation (mask)
-                                                                                // so renormalizationFactor is not initialized
-  }
-}
-
-}
+  if (detectMissingValue)
+  {
+    double defaultValue = std::numeric_limits<double>::quiet_NaN();
+    size_t nbIndex=dataOut.numElements() ; 
+
+    for (int idx = 0; idx < nbIndex; ++idx)
+    {
+      if (allMissing(idx)) dataOut(idx) = defaultValue; // If all data source are nan then data destination must be nan
+    }
+    
+    if (renormalize)
+    {
+      if (renormalizationFactor.numElements()>0) dataOut/=renormalizationFactor ; // In some case, process doesn't received any data for interpolation (mask)
+                                                                                 // so renormalizationFactor is not initialized
+    }
+  }
+}
+
+}

XIOS/dev/XIOS_DEV_CMIP6/src/transformation/domain_algorithm_interpolate.hpp

@@ -55 +55 @@
 private:
   CArray<double,1> renormalizationFactor ;
+  CArray<bool,1> allMissing ;
   bool detectMissingValue ;
   bool renormalize ;

XIOS/dev/XIOS_DEV_CMIP6/src/transformation/generic_algorithm_transformation.cpp

@@ -37 +37 @@
   int nbLocalIndex = localIndex.size();   
   double defaultValue = std::numeric_limits<double>::quiet_NaN();
+    
   if (ignoreMissingValue)
   {
+    if (firstPass) dataOut=defaultValue ;
+    
     for (int idx = 0; idx < nbLocalIndex; ++idx)
     {
-      if (NumTraits<double>::isNan(*(dataInput + idx)))
-      {
-        flagInitial[localIndex[idx].first] = false;
-      }
-      else
-      {
-        dataOut(localIndex[idx].first) += *(dataInput + idx) * localIndex[idx].second;
-        flagInitial[localIndex[idx].first] = true; // Reset flag to indicate not all data source are nan
-      }
-    }
-
-    // If all data source are nan then data destination must be nan
-    for (int idx = 0; idx < nbLocalIndex; ++idx)
-    {
-      if (!flagInitial[localIndex[idx].first])
-        dataOut(localIndex[idx].first) = defaultValue;
-    }
+      if (! NumTraits<double>::isNan(*(dataInput + idx)))
+      {
+        if (flagInitial[localIndex[idx].first]) dataOut(localIndex[idx].first) = *(dataInput + idx) * localIndex[idx].second;
+        else dataOut(localIndex[idx].first) += *(dataInput + idx) * localIndex[idx].second;
+        flagInitial[localIndex[idx].first] = false; // Reset flag to indicate not all data source are nan
+      }
+    }
+
   }
   else