Merge pull request #166 from KNMI/fix-postproc-hdf5-data

Fix postproc hdf5 data

CCDFDataModel/CCDFVariable.cpp

@@ -127,16 +127,15 @@ int CDF::Variable::readData(CDFType readType,size_t *_start,size_t *_count,ptrdi
       float fscale = float(scaleFactor);
       float foffset = float(addOffset);
       for(size_t j=0;j<lsize;j++)scaleData[j]=scaleData[j]*fscale+foffset;
-      fillValue=fillValue*fscale+foffset;
-      float f=(float)fillValue;
-      if( hasFillValue)getAttribute("_FillValue")->setData(CDF_FLOAT,&f,1);
+      float newFillValue = fillValue * fscale + foffset;
+      if( hasFillValue)getAttribute("_FillValue")->setData(CDF_FLOAT,&newFillValue,1);
     }
 
     if(scaleType == CDF_DOUBLE){
       double *scaleData = (double*)data;
       for(size_t j=0;j<lsize;j++)scaleData[j]=scaleData[j]*scaleFactor+addOffset;
-      fillValue=fillValue*scaleFactor+addOffset;
-      if( hasFillValue)getAttribute("_FillValue")->setData(CDF_DOUBLE,&fillValue,1);
+      double newFillValue=fillValue*scaleFactor+addOffset;
+      if( hasFillValue)getAttribute("_FillValue")->setData(CDF_DOUBLE,&newFillValue,1);
     }
     //removeAttribute("scale_factor");
     //removeAttribute("add_offset");

Dockerfile

@@ -4,7 +4,7 @@ USER root
 LABEL maintainer="[email protected]"
 
 # Version should be same as in Definitions.h
-LABEL version="2.5.11" 
+LABEL version="2.5.12" 
 
 ######### First stage (build) ############
 

adagucserverEC/CDataPostProcessor.cpp

@@ -791,6 +791,50 @@ int CDPPExecutor::executeProcessors( CDataSource *dataSource,int mode){
   return 0;
 }
 
+int CDPPExecutor::executeProcessors( CDataSource *dataSource,int mode, double * data, size_t numItems){
+  //const CT::PointerList<CDPPInterface*> *availableProcs = getPossibleProcessors();
+  //CDBDebug("executeProcessors, found %d",dataSource->cfgLayer->DataPostProc.size());
+  for(size_t dpi=0;dpi<dataSource->cfgLayer->DataPostProc.size();dpi++){
+    CServerConfig::XMLE_DataPostProc * proc = dataSource->cfgLayer->DataPostProc[dpi];
+    for(size_t procId = 0;procId<dataPostProcessorList->size();procId++){
+      int code = dataPostProcessorList->get(procId)->isApplicable(proc,dataSource);
+
+      if(code == CDATAPOSTPROCESSOR_CONSTRAINTSNOTMET){
+        CDBError("Constraints for DPP %s are not met",dataPostProcessorList->get(procId)->getId());
+      }
+
+      /*Will be runned when datasource metadata been loaded */
+      if(mode == CDATAPOSTPROCESSOR_RUNBEFOREREADING){
+        if(code&CDATAPOSTPROCESSOR_RUNBEFOREREADING){
+          try{
+            int status = dataPostProcessorList->get(procId)->execute(proc,dataSource,CDATAPOSTPROCESSOR_RUNBEFOREREADING,  NULL, 0);
+            if(status != 0){
+              CDBError("Processor %s failed RUNBEFOREREADING, statuscode %d",dataPostProcessorList->get(procId)->getId(),status);
+            }
+          }catch(int e){
+            CDBError("Exception in Processor %s failed RUNBEFOREREADING, exceptioncode %d",dataPostProcessorList->get(procId)->getId(),e);
+          }
+        }
+      }
+
+      /*Will be runned when datasource data been loaded */
+      if(mode == CDATAPOSTPROCESSOR_RUNAFTERREADING){
+        if(code&CDATAPOSTPROCESSOR_RUNAFTERREADING){
+          try{
+            int status = dataPostProcessorList->get(procId)->execute(proc,dataSource,CDATAPOSTPROCESSOR_RUNAFTERREADING,  data, numItems);
+            if(status != 0){
+              CDBError("Processor %s failed RUNAFTERREADING, statuscode %d",dataPostProcessorList->get(procId)->getId(),status);
+            }
+          }catch(int e){
+            CDBError("Exception in Processor %s failed RUNAFTERREADING, exceptioncode %d",dataPostProcessorList->get(procId)->getId(),e);
+          }
+        }
+      }
+    }
+  }
+  return 0;
+}
+
 
 /************************/
 /*      CDPPBeaufort     */
@@ -1114,15 +1158,36 @@ int CDPDBZtoRR::isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource
 }
 
 float CDPDBZtoRR::getRR(float dbZ) {
+  // TODO: Check why -32 as input does not return 0
   return pow((pow(10,dbZ/10.)/200),1/1.6);
 }
 
+int CDPDBZtoRR::execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){
+  CDBDebug("CDPDBZtoRR");
+  if((isApplicable(proc,dataSource)&mode)==false){
+    return -1;
+  }
+  if(mode==CDATAPOSTPROCESSOR_RUNBEFOREREADING){  
+    dataSource->getDataObject(0)->setUnits("mm/hr");
+  }
+  if(mode==CDATAPOSTPROCESSOR_RUNAFTERREADING){  
+    float noDataValue=dataSource->getDataObject(0)->dfNodataValue;
+    for(size_t j = 0; j < numItems; j++) {
+      if (data[j]==data[j]) { // Check if NaN
+        if (data[j]!=noDataValue) { // Check if equal to nodatavalue of the datasource
+          data[j]=getRR(data[j]);
+        }
+      }
+    }
+  }
+  return 0;
+}
+
 int CDPDBZtoRR::execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode){
   if((isApplicable(proc,dataSource)&mode)==false){
     return -1;
   }
   if(mode==CDATAPOSTPROCESSOR_RUNBEFOREREADING){  
-//    dataSource->getDataObject(0)->cdfVariable->setAttributeText("units","mm/hr");
     dataSource->getDataObject(0)->setUnits("mm/hr");
   }
   if(mode==CDATAPOSTPROCESSOR_RUNAFTERREADING){  

adagucserverEC/CDataPostProcessor.h

@@ -51,9 +51,14 @@
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource) = 0;
 
     /**
-    * Executes the data postprocessor
+    * Executes the data postprocessor on a datasource, on the full grid
     */
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode) = 0;
+
+    /**
+    * Executes the data postprocessor for a given array
+    */
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems) = 0;
   };
 
   /**
@@ -66,6 +71,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
   /**
@@ -94,6 +100,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
 
@@ -114,6 +121,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
   /**
@@ -145,6 +153,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
   /**
@@ -181,6 +190,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
   class CDPPExecutor{
@@ -192,6 +202,7 @@
     ~CDPPExecutor();
     const CT::PointerList<CDPPInterface*> *getPossibleProcessors();
     int executeProcessors( CDataSource *dataSource,int mode);
+    int executeProcessors( CDataSource *dataSource,int mode, double * data, size_t numItems);
   };
 
   /**
@@ -205,6 +216,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
   /**
@@ -217,6 +229,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
   /**
@@ -230,6 +243,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems);
   };
 
   /**
@@ -243,6 +257,7 @@
     virtual const char *getId();
     virtual int isApplicable(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource);
     virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode);
+    virtual int execute(CServerConfig::XMLE_DataPostProc* proc, CDataSource* dataSource,int mode, double *data, size_t numItems){return 1;} // TODO: Still need to implement
   };
 
 

adagucserverEC/CMakeJSONTimeSeries.cpp

@@ -461,7 +461,6 @@ class UniqueRequests{
         if(projCacheInfo.isOutsideBBOX == false){
           CDFObject *cdfObject = CDFObjectStore::getCDFObjectStore()->getCDFObjectHeader(dataSource, dataSource->srvParams,(filemapiterator->first).c_str());
 
-
           if(cdfObject->getVariableNE("forecast_reference_time")!=NULL){
             CDBDebug("IS REFERENCE TIME");
             CDF::Dimension *forecastRefDim = new CDF::Dimension();
@@ -552,11 +551,19 @@ class UniqueRequests{
               CDBDebug("  %d %s [%d:%d]",i,"",start[i],count[i]);
             }
 #endif       
-
-  /*
-            if(variable->getAttributeNE("scale_factor")!=NULL){
-              readDataAsCDFDouble = true;
-            }*/
+
+
+            /* 
+             * In case a scale_factor and add_offset attribute is present, we need to read the data into the same datatype as this attribute
+             * This allows it to be unpacked properly to the final scaled values 
+             */
+            CDF::Attribute *scale_factor = variable->getAttributeNE("scale_factor");
+            if(scale_factor!=NULL) {
+              variable->setType(CDF_FLOAT);
+              if(scale_factor->getType()==CDF_DOUBLE){
+                variable->setType(CDF_DOUBLE);
+              }
+            }
 
             if(readDataAsCDFDouble){
               variable->setType(CDF_DOUBLE);
@@ -572,7 +579,7 @@ class UniqueRequests{
               /**
               * DataPostProc: Here our datapostprocessor comes into action!
               */
-              for(size_t dpi=0;dpi<dataSource->cfgLayer->DataPostProc.size();dpi++){
+              for(size_t dpi=0;dpi<dataSource->cfgLayer->DataPostProc.size();dpi++) {
                 CServerConfig::XMLE_DataPostProc * proc = dataSource->cfgLayer->DataPostProc[dpi];
                 //Algorithm ax+b:
                 if(proc->attr.algorithm.equals("ax+b")){
@@ -594,8 +601,10 @@ class UniqueRequests{
                 //Apply units:
                 if(proc->attr.units.empty()==false){
                   dataSource->getDataObject(dataObjectNr)->setUnits(proc->attr.units.c_str());
-                }
-
+                } 
+              }
+              if (readDataAsCDFDouble) {
+                CDataPostProcessor::getCDPPExecutor()->executeProcessors(dataSource,CDATAPOSTPROCESSOR_RUNAFTERREADING, (double*)variable->data, variable->getSize());
               }
               /* End of data postproc */
 #ifdef CMakeJSONTimeSeries_DEBUG          
@@ -659,14 +668,17 @@ int CMakeJSONTimeSeries::MakeJSONTimeSeries(CDrawImage *drawImage,CImageWarper *
   /**
   * DataPostProc: Here our datapostprocessor comes into action!
   */
-  for(size_t dpi=0;dpi<dataSource->cfgLayer->DataPostProc.size();dpi++){
-    CServerConfig::XMLE_DataPostProc * proc = dataSource->cfgLayer->DataPostProc[dpi];
-    //Algorithm ax+b:
-    if(proc->attr.algorithm.equals("ax+b")){
-      uniqueRequest.readDataAsCDFDouble = true;
-      break;
-    }
+  if (dataSource->cfgLayer->DataPostProc.size() > 0) {
+    uniqueRequest.readDataAsCDFDouble = true;
   }
+  // for(size_t dpi=0;dpi<dataSource->cfgLayer->DataPostProc.size();dpi++){
+  //   CServerConfig::XMLE_DataPostProc * proc = dataSource->cfgLayer->DataPostProc[dpi];
+  //   //Algorithm ax+b:
+  //   if(proc->attr.algorithm.equals("ax+b")){
+  //     uniqueRequest.readDataAsCDFDouble = true;
+  //     break;
+  //   }
+  // }
 
   int numberOfDims = dataSource->requiredDims.size();
   int numberOfSteps = dataSource->getNumTimeSteps();

adagucserverEC/Definitions.h

@@ -28,7 +28,7 @@
 #ifndef Definitions_H
 #define Definitions_H
 
-#define ADAGUCSERVER_VERSION "2.5.11" // Please also update in the Dockerfile to the same version
+#define ADAGUCSERVER_VERSION "2.5.12" // Please also update in the Dockerfile to the same version
 
 //CConfigReaderLayerType
 #define CConfigReaderLayerTypeUnknown 0

data/config/datasets/adaguc.KNMIHDF5.test.xml

@@ -0,0 +1,51 @@
+<?xml version="1.0" encoding="UTF-8" ?>
+<Configuration>
+
+  <WMS>
+    <Title>RADNL_OPER_R___25PCPRR_L3_WMS</Title>
+    <Abstract>Radar precipitation measurements above the Netherlands, on a 1.0x1.0 km grid, measurements are available in a five minute time interval. The intensity is in kg/m2/hour (mm/hour). The dataset is created from KNMI RAD_NL25_PCP_NA files. For interoperability, the original unit reflectivity in DBZ is converted to precipitation flux in kg/m2/h. The conversion from dbZ to kg/m2/h is applied with the formula R = 10^((PixelValue -109)/32).</Abstract>
+    <RootLayer>
+      <Title>RADNL_OPER_R___25PCPRR_L3_WMS</Title>
+      <Abstract></Abstract>
+    </RootLayer>
+  </WMS>
+
+  <WCS>
+    <Title>RADNL_OPER_R___25PCPRR_L3_WCS</Title>
+    <Label>wcsLabel</Label>
+    <WCSFormat name="netcdf" driver="ADAGUCNetCDF" mimetype="Content-Type:application/netcdf" options=""/>
+    <WCSFormat name="NetCDF3" title="GDAL NetCDF3" driver="ADAGUC"   mimetype="Content-Type:text/plain" options="FORCENC3=TRUE"/>
+    <WCSFormat name="NetCDF4" driver="ADAGUC"   mimetype="Content-Type:text/plain"/>
+  </WCS>
+
+
+  <Legend name="KNMIRadarPalette_mmh" type="interval">
+    <palette min="0"  max="38"  red="255" green="255" blue="255"/>
+    <palette min="38"  max="80" red="170" green="170" blue="170"/>
+    <palette min="80" max="118" red="85"  green="85"  blue="85"/>
+    <palette min="118" max="160" red="255" green="128" blue="128"/>
+    <palette min="160" max="198" red="255" green="0" blue="0"/>
+    <palette min="198" max="255" red="0"   green="0"   blue="0"/>
+  </Legend>
+
+ <Style name="radar">
+    <Legend>KNMIRadarPalette_mmh</Legend>
+    <ValueRange min=".09" max="1000000"/>
+    <RenderMethod>nearest</RenderMethod>
+    <Min>0.1</Min>
+    <Max>100</Max>
+    <Log>10</Log>
+  </Style>
+
+  <!-- Layers -->
+  <Layer type="database" hidden="false">
+    <Name>RAD_NL25_PCP_CM</Name>
+    <Title>Precipitation Radar NL</Title>
+    <Variable>image1.image_data</Variable>
+    <!--               RAD_NL25_PCP_CM_202102170000.h5 -->
+    <FilePath filter="^RAD_NL25_PCP_CM_.*\.h5$">{ADAGUC_PATH}/data/datasets/test/</FilePath>
+    <Dimension name="time" units="ISO8601" interval="PT5M" default="max">time</Dimension>
+    <Styles>radar</Styles>
+    <DataPostProc algorithm="dbztorr" units="mm/hr"/>
+  </Layer>
+</Configuration>

tests/AdagucTests/TestWMS.py

@@ -574,3 +574,17 @@ def test_WMSGetMapCustomCRSEPSG3413Projection_ipcc_cmip5_tas_historical_subset_n
     #     self.assertEqual(status, 0)
     #     self.assertEqual(data.getvalue(), AdagucTestTools().readfromfile(self.expectedoutputsspath + filename))
 
+    def test_WMSGetFeatureInfo_timeseries_KNMIHDF5_json(self):
+        AdagucTestTools().cleanTempDir()
+        ADAGUC_PATH = os.environ['ADAGUC_PATH']
+        env={'ADAGUC_CONFIG' : ADAGUC_PATH + "/data/config/adaguc.tests.dataset.xml," + ADAGUC_PATH + "/data/config/datasets/adaguc.KNMIHDF5.test.xml"}
+        config =  ADAGUC_PATH + '/data/config/adaguc.tests.dataset.xml,' + ADAGUC_PATH + '/data/config/datasets/adaguc.KNMIHDF5.test.xml'
+        status,data,headers = AdagucTestTools().runADAGUCServer(args = ['--updatedb', '--config', config], env = self.env, isCGI = False)
+        self.assertEqual(status, 0)
+
+        filename="test_WMSGetFeatureInfo_timeseries_KNMIHDF5_json.json"
+        status,data,headers = AdagucTestTools().runADAGUCServer("dataset=adaguc.KNMIHDF5.test&service=WMS&request=GetFeatureInfo&version=1.3.0&layers=RAD_NL25_PCP_CM&query_layers=RAD_NL25_PCP_CM&crs=EPSG%3A3857&bbox=467411.5837657447%2C5796421.971094566%2C889884.3758374067%2C7834481.671540775&width=199&height=960&i=103&j=501&format=image%2Fgif&info_format=application%2Fjson&time=1000-01-01T00%3A00%3A00Z%2F3000-01-01T00%3A00%3A00Z&", env = env)
+        AdagucTestTools().writetofile(self.testresultspath + filename,data.getvalue())
+        self.assertEqual(status, 0)
+        self.assertEqual(data.getvalue(), AdagucTestTools().readfromfile(self.expectedoutputsspath + filename))
+

tests/expectedoutputs/TestWMS/test_WMSGetFeatureInfo_timeseries_KNMIHDF5_json.json

@@ -0,0 +1 @@
+[{"name":"RAD_NL25_PCP_CM","standard_name":"image1.image_data","units":"mm/hr","point":{"SRS":"EPSG:4326","coords":"6.163147,51.837904"},"dims":"time","data":{"2021-06-22T20:00:00Z":"5.225240"}}]