Merge commit '631749fc0d93de60626c1ea67240068613cec866'

agrolib/dbMeteoGrid/dbMeteoGrid.cpp

@@ -4126,7 +4126,7 @@ bool Crit3DMeteoGridDbHandler::MeteoGridToRasterFlt(double cellSize, const gis::
             {
                 myGrid.getXY(row, col, utmx, utmy);
                 gis::getLatLonFromUtm(gisSettings, utmx, utmy, &lat, &lon);
-                gis::getGridRowColFromXY (latlonHeader, lon, lat, &dataGridRow, &dataGridCol);
+                gis::getGridRowColFromLonLat (latlonHeader, lon, lat, &dataGridRow, &dataGridCol);
                 if (dataGridRow < 0 || dataGridRow >= latlonHeader.nrRows || dataGridCol < 0 || dataGridCol >= latlonHeader.nrCols)
                 {
                     myValue = NODATA;

agrolib/gis/color.cpp

@@ -56,6 +56,7 @@ Crit3DColorScale::Crit3DColorScale()
     _maximum = NODATA;
     _isFixedRange = false;
     _isHideOutliers = false;
+    _isTransparent = false;
 
     _classification = classificationMethod::EqualInterval;
 }
@@ -212,8 +213,8 @@ bool setSlopeStabilityScale(Crit3DColorScale* myScale)
     myScale->keyColor[0] = Crit3DColor(128, 0, 128);       /*!< violet */
     myScale->keyColor[1] = Crit3DColor(255, 0, 0);         /*!< red */
     myScale->keyColor[2] = Crit3DColor(255, 255, 0);       /*!< yellow */
-    myScale->keyColor[3] = Crit3DColor(196, 196, 32);      /*!< yellow/green */
-     myScale->keyColor[4] = Crit3DColor(64, 196, 64);      /*!< green */
+    myScale->keyColor[3] = Crit3DColor(64, 196, 64);       /*!< green */
+    myScale->keyColor[4] = Crit3DColor(128, 255, 128);     /*!< light green */
 
     return(myScale->classify());
 }
@@ -365,6 +366,23 @@ bool reverseColorScale(Crit3DColorScale* myScale)
 }
 
 
+void mixColor(const Crit3DColor &backColor, const Crit3DColor &foreColor, Crit3DColor &colorOut, float alpha)
+{
+    if (alpha == 0)
+    {
+        colorOut.red = backColor.red;
+        colorOut.green = backColor.green;
+        colorOut.blue = backColor.blue;
+    }
+    else
+    {
+        colorOut.red = std::min(255, int(backColor.red * (1. - alpha)) + int(foreColor.red * alpha));
+        colorOut.green = std::min(255, int(backColor.green * (1. - alpha)) + int(foreColor.green * alpha));
+        colorOut.blue = std::min(255, int(backColor.blue * (1. - alpha)) +int( foreColor.blue * alpha));
+    }
+}
+
+
 /*!
  * \brief roundColorScale round colorScale values on the second (or third) digit of each range.
  * It requires that nrColors is a multiply of nrIntervals for a correct visualization in the colors legend.

agrolib/gis/color.h

@@ -15,7 +15,7 @@
         short blue;
 
         Crit3DColor();
-        Crit3DColor(short,short,short);
+        Crit3DColor(short, short, short);
     };
 
     class Crit3DColorScale {
@@ -26,6 +26,7 @@
         double _minimum, _maximum;
         bool _isFixedRange;
         bool _isHideOutliers;
+        bool _isTransparent;
         int _classification;
 
     public:
@@ -55,6 +56,9 @@
 
         void setHideOutliers(bool hideOutliers) { _isHideOutliers = hideOutliers; }
         bool isHideOutliers() { return _isHideOutliers; }
+
+        void setTransparent(bool transparent) { _isTransparent = transparent; }
+        bool isTransparent() { return _isTransparent; }
     };
 
     bool setDefaultScale(Crit3DColorScale* myScale);
@@ -75,5 +79,7 @@
     bool setSurfaceWaterScale(Crit3DColorScale* myScale);
     bool setLAIScale(Crit3DColorScale* myScale);
 
+    void mixColor(const Crit3DColor &backColor, const Crit3DColor &foreColor, Crit3DColor &colorOut, float alpha);
+
 
 #endif // CRIT3DCOLOR_H

agrolib/gis/gis.cpp

@@ -791,7 +791,7 @@ namespace gis
         v->col = int(floor((p.x - myHeader.llCorner.x) / myHeader.cellSize));
     }
 
-    void getGridRowColFromXY (const Crit3DLatLonHeader& myHeader, double myX, double myY, int *row, int *col)
+    void getGridRowColFromLonLat(const Crit3DLatLonHeader& myHeader, double myX, double myY, int *row, int *col)
     {
         *row = int(floor((myY - myHeader.llCorner.latitude) / myHeader.dy));
         *col = int(floor((myX - myHeader.llCorner.longitude) / myHeader.dx));

agrolib/gis/gis.h

@@ -217,7 +217,7 @@
         void getRowColFromXY(const Crit3DRasterHeader& myHeader, double myX, double myY, int *row, int *col);
         void getRowColFromXY(const Crit3DRasterHeader& myHeader, const Crit3DUtmPoint& p, int *row, int *col);
         void getRowColFromXY(const Crit3DRasterHeader& myHeader, const Crit3DUtmPoint& p, Crit3DRasterCell* v);
-        void getGridRowColFromXY(const Crit3DLatLonHeader& myHeader, double myX, double myY, int *row, int *col);
+        void getGridRowColFromLonLat(const Crit3DLatLonHeader& myHeader, double myX, double myY, int *row, int *col);
 
         void getRowColFromLatLon(const Crit3DLatLonHeader &latLonHeader, const Crit3DGeoPoint& p, int *myRow, int *myCol);
         bool isOutOfGridRowCol(int myRow, int myCol, const Crit3DRasterGrid &rasterGrid);

agrolib/grapevine/grapevine.cpp

@@ -113,11 +113,6 @@ void Vine3D_Grapevine::resetLayers()
 {
     for (int i=0 ; i < nrMaxLayers ; i++)
     {
-        //psiSoilProfile[i] = NODATA ;
-        //soilWaterContentProfile[i]= NODATA ;
-        //soilWaterContentProfileFC[i]= NODATA;
-        //soilWaterContentProfileWP[i]= NODATA;
-        //soilFieldCapacity[i] = NODATA;
         fractionTranspirableSoilWaterProfile[i]= NODATA;
         stressCoefficientProfile[i] = NODATA;
         transpirationInstantLayer[i] = NODATA;
@@ -126,15 +121,10 @@ void Vine3D_Grapevine::resetLayers()
     }
 }
 
-bool Vine3D_Grapevine::initializeLayers(int myMaxLayers)
+void Vine3D_Grapevine::initializeLayers(int myMaxLayers)
 {
     nrMaxLayers = myMaxLayers;
 
-    //psiSoilProfile = (double *) calloc(nrLayers, sizeof(double));
-    //soilWaterContentProfile = (double *) calloc(nrLayers, sizeof(double));
-    //soilWaterContentProfileFC = (double *) calloc(nrLayers, sizeof(double));
-    //soilWaterContentProfileWP = (double *) calloc(nrLayers, sizeof(double));
-    //soilFieldCapacity = (double *) calloc (nrLayers, sizeof(double));
     fractionTranspirableSoilWaterProfile = static_cast<double*> (calloc(size_t(nrMaxLayers), sizeof(double)));
     stressCoefficientProfile = static_cast<double*> (calloc(size_t(nrMaxLayers), sizeof(double)));
     transpirationInstantLayer = static_cast<double*> (calloc(size_t(nrMaxLayers), sizeof(double)));
@@ -143,10 +133,9 @@ bool Vine3D_Grapevine::initializeLayers(int myMaxLayers)
     currentProfile = static_cast<double*> (calloc(size_t(nrMaxLayers), sizeof(double)));
 
     resetLayers();
-
-    return true;
 }
 
+
 void Vine3D_Grapevine::setDate (Crit3DTime myTime)
 {
     myDoy = getDoyFromDate(myTime.date);

agrolib/grapevine/grapevine.h

@@ -192,13 +192,13 @@
     struct Crit3DModelCase {
         int id;
         Crit3DLanduse landuse;
-        //int soilIndex;
+        int soilIndex;
 
         float shootsPerPlant;
         float plantDensity;
         float maxLAIGrass;
         int trainingSystem;
-        float maxIrrigationRate;        //[mm/h]
+        float maxIrrigationRate;        // [mm h-1]
 
         int soilLayersNr;
         double soilTotalDepth;
@@ -393,7 +393,7 @@
         Vine3D_Grapevine();
 
         //void initializeGrapevineModel(TVineCultivar* cultivar, double secondsPerStep);
-        bool initializeLayers(int myMaxLayers);
+        void initializeLayers(int myMaxLayers);
         bool initializeStatePlant(int doy, Crit3DModelCase *vineField);
         void resetLayers();
 

agrolib/graphics/mapGraphicsRasterObject.cpp

@@ -519,7 +519,7 @@ bool RasterObject::getRowCol(gis::Crit3DGeoPoint geoPoint, int* row, int* col)
     if (! this->isGrid)
         return false;
 
-    gis::getGridRowColFromXY(this->latLonHeader, geoPoint.longitude, geoPoint.latitude, row, col);
+    gis::getGridRowColFromLonLat(this->latLonHeader, geoPoint.longitude, geoPoint.latitude, row, col);
 
     // check out of grid
     if (gis::isOutOfGridRowCol(*row, *col, this->latLonHeader))

agrolib/interpolation/interpolation.cpp

@@ -1068,10 +1068,11 @@ void localSelection(vector <Crit3DInterpolationDataPoint> &inputPoints, vector <
         for (i=0; i< selectedPoints.size(); i++)
         {
             //selectedPoints[i].regressionWeight = MAXVALUE(1 - selectedPoints[i].distance / (std::pow(maxDistance, 16.0/15.0)), EPSILON);
-            selectedPoints[i].regressionWeight = MAXVALUE(std::exp(-std::pow(selectedPoints[i].distance/(0.8*maxDistance),3.0)),EPSILON);
+            //selectedPoints[i].regressionWeight = MAXVALUE(std::exp(-std::pow(selectedPoints[i].distance/((4/5*maxDistance)),7.0)),EPSILON);
+            selectedPoints[i].regressionWeight = MAXVALUE((-(1/std::pow(maxDistance,4)*(std::pow(selectedPoints[i].distance,4)))+1),EPSILON);
             //selectedPoints[i].regressionWeight = 1;
-            selectedPoints[i].heightWeight = 1./((2./maxHeightDelta)*selectedPoints[i].point->z+1);
-            //selectedPoints[i].heightWeight = 1;
+            //selectedPoints[i].heightWeight = 1./((2./maxHeightDelta)*selectedPoints[i].point->z+1);
+            selectedPoints[i].heightWeight = 1;
         }
     mySettings.setLocalRadius(float(maxDistance));
 }
@@ -1689,9 +1690,9 @@ bool multipleDetrendingMain(std::vector <Crit3DInterpolationDataPoint> &myPoints
         elevationCombination.resetCombination(mySettings->getSelectedCombination().getProxySize());
         elevationCombination.setProxyActive(elevationPos, true);
 
-        if (parameters.empty()) {
+        if (parameters.empty())
             parameters.resize(elevationPos + 1);
-        }
+
 
         if (!multipleDetrendingElevation(elevationCombination, parameters[elevationPos], myPoints, mySettings, myVar, errorStr))
             return false;
@@ -1871,7 +1872,7 @@ bool multipleDetrendingElevation(Crit3DProxyCombination elevationCombination, st
     interpolation::bestFittingMarquardt_nDimension_singleFunction(*func, nrMaxStep, 4, parametersMin[elevationPos], parametersMax[elevationPos], parameters[elevationPos], parametersDelta[elevationPos],
                                                    100, 0.005, 0.002, predictors, predictands, weights);
 
-    mySettings->setSingleFittingFunction(myFunc[elevationPos], elevationPos);
+
     mySettings->setSingleFittingParameters(parameters[elevationPos], elevationPos);
 
 
@@ -2104,7 +2105,7 @@ bool multipleDetrending(Crit3DProxyCombination othersCombination, std::vector<st
     //std::vector <std::vector<double>> parameters;
     std::vector<std::function<double(double, std::vector<double>&)>> myFunc(mySettings->getCurrentCombination().getProxySize(), nullptr);
 
-    unsigned int nrMaxStep = 100;
+    unsigned int nrMaxStep = 10;
     if (parameters.empty())
         nrMaxStep *= 10;
 
@@ -2130,7 +2131,7 @@ bool multipleDetrending(Crit3DProxyCombination othersCombination, std::vector<st
     interpolation::bestFittingMarquardt_nDimension(&functionSum, myFunc, nrMaxStep, 4, parametersMin, parametersMax, parameters, parametersDelta,
                                                    100, 0.005, 0.002, predictors, predictands, weights, elevationPos);
 
-
+    myFunc.clear();
     int k = 0;
     for (unsigned int i = 0; i < othersCombination.getProxySize(); i++)
     {
@@ -2142,11 +2143,10 @@ bool multipleDetrending(Crit3DProxyCombination othersCombination, std::vector<st
                 mySettings->setSingleFittingFunction(functionLinear, i);
             }
             k++;
-            myFunc[i] = functionLinear;
+            myFunc.push_back(functionLinear);
         }
     }
 
-
     std::vector <double> proxyValues;
 
     // detrending

agrolib/interpolation/interpolationSettings.cpp

@@ -301,7 +301,7 @@ void Crit3DInterpolationSettings::setFittingParameters(const std::vector<std::ve
 void Crit3DInterpolationSettings::setSingleFittingParameters(std::vector<double> &newFittingParameters, int paramIndex)
 {
     if (fittingParameters.size() <= paramIndex)
-            fittingParameters.resize(paramIndex+1);
+        fittingParameters.resize(paramIndex+1);
     fittingParameters[paramIndex] = newFittingParameters;
 }
 

agrolib/mathFunctions/furtherMathFunctions.cpp

@@ -1402,27 +1402,14 @@ namespace interpolation
                             }
                             counter++;
 
-                            if (!parameters.size() > proxyIndex && !parameters[proxyIndex].empty())
+                            if (parameters.size() > proxyIndex && !parameters[proxyIndex].empty())
                             {
                                 if (dir == 0)
                                     parameters[proxyIndex][paramIndex] = MINVALUE(firstGuessParam[proxyIndex][paramIndex] + directions[dir] * step * steps[paramIndex], parametersMax[proxyIndex][paramIndex]);
                                 else
                                     parameters[proxyIndex][paramIndex] = MAXVALUE(firstGuessParam[proxyIndex][paramIndex] + directions[dir] * step * steps[paramIndex], parametersMin[proxyIndex][paramIndex]);
                             }
 
-                            /*for (i=0; i<nrPredictors; i++)
-                    {
-                        for (j=0; j<nrParameters[i]; j++)
-                        {
-                            do {
-                                truncNormal = normal_dis(gen);
-                            } while(truncNormal <= 0.0 || truncNormal >= 1.0);
-                            parameters[i][j] = parametersMin[i][j] + (truncNormal)*(parametersMax[i][j]-parametersMin[i][j]);
-                        }
-                    }
-*/
-
-
                             if ((counter > nrTrials) || ((R2Previous[0] != NODATA) && fabs(R2Previous[0]-R2Previous[nrMinima-1]) < deltaR2 ))
                             {
                                 for (i=0;i<parameters.size();i++)
@@ -2082,7 +2069,7 @@ namespace interpolation
 
         //grigliato
         std::vector<double> steps;
-        const int numSteps = 50;
+        const int numSteps = 30;
         if (parameters.size() == 4)
             steps = {2*(parametersMax[0]-parametersMin[0])/numSteps, 2*(parametersMax[1]-parametersMin[1])/numSteps, 2*(parametersMax[2]-parametersMin[2])/numSteps, 2*(parametersMax[3]-parametersMin[3])/numSteps};
         else if (parameters.size() == 6)

agrolib/meteo/meteo.cpp

@@ -790,6 +790,7 @@ bool setColorScale(meteoVariable variable, Crit3DColorScale *colorScale)
                 || variable == snowLiquidWaterContent || variable == snowMelt)
             {
                 colorScale->setHideOutliers(true);
+                colorScale->setTransparent(true);
             }
             break;  
         case snowAge:

agrolib/meteo/meteoGrid.cpp

@@ -988,7 +988,7 @@ void Crit3DMeteoGrid::saveRowColfromZone(gis::Crit3DRasterGrid* zoneGrid, std::v
                     double utmX = x;
                     double utmY = y;
                     gis::getLatLonFromUtm(_gisSettings, utmX, utmY, &y, &x);
-                    gis::getGridRowColFromXY(_gridStructure.header(), x, y, &myRow, &myCol);
+                    gis::getGridRowColFromLonLat(_gridStructure.header(), x, y, &myRow, &myCol);
                 }
                 else
                 {

agrolib/project/project.cpp

@@ -3230,7 +3230,7 @@ bool Project::loadProject()
     if (dbPointsFileName != "")
         if (! loadMeteoPointsDB(dbPointsFileName))
         {
-            errorString = "load Meteo Points DB failed";
+            errorString = "load Meteo Points DB failed:\n" + dbPointsFileName;
             errorType = ERROR_DBPOINT;
             logError();
             return false;
@@ -3709,7 +3709,7 @@ void Project::showLocalProxyGraph(gis::Crit3DGeoPoint myPoint, gis::Crit3DRaster
     }
     if (this->meteoGridLoaded && !this->meteoGridDbHandler->meteoGrid()->dataMeteoGrid.parametersCell.empty())
     {
-        gis::getGridRowColFromXY(meteoGridDbHandler->meteoGrid()->gridStructure().header(), myPoint.longitude, myPoint.latitude, &row, &col);
+        gis::getGridRowColFromLonLat(meteoGridDbHandler->meteoGrid()->gridStructure().header(), myPoint.longitude, myPoint.latitude, &row, &col);
         parameters = meteoGridDbHandler->meteoGrid()->dataMeteoGrid.getParametersFromRowCol(row, col);
     }
 
@@ -4200,11 +4200,16 @@ bool Project::exportMeteoGridToCsv(QString fileName)
         for (int col = 0; col < meteoGridDbHandler->meteoGrid()->dataMeteoGrid.header->nrCols; col++)
         {
             float value = meteoGridDbHandler->meteoGrid()->dataMeteoGrid.value[row][col];
-            std::string id = meteoGridDbHandler->meteoGrid()->meteoPoints()[row][col]->id;
-            std::string name = meteoGridDbHandler->meteoGrid()->meteoPoints()[row][col]->name;
 
             if (value != NO_ACTIVE && value != NODATA)
             {
+                int newRow = row;
+                if (! meteoGridDbHandler->meteoGrid()->gridStructure().isUTM())
+                    newRow = meteoGridDbHandler->meteoGrid()->gridStructure().nrRow() - 1 - row;
+
+                std::string id = meteoGridDbHandler->meteoGrid()->meteoPoints()[newRow][col]->id;
+                std::string name = meteoGridDbHandler->meteoGrid()->meteoPoints()[newRow][col]->name;
+
                 out << QString::fromStdString(id + ',' + name + ',') + QString::number(value) + "\n";
             }
         }

agrolib/soilFluxes3D/soilFluxes3D.cpp

@@ -658,7 +658,7 @@ int DLL_EXPORT __STDCALL setHydraulicProperties(int waterRetentionCurve,
 
  /*!
   * \brief getTotalWaterContent
-  * \return total water content [m^3]
+  * \return total water content [m3]
   */
  double DLL_EXPORT __STDCALL getTotalWaterContent()
  {
@@ -823,10 +823,15 @@ int DLL_EXPORT __STDCALL setHydraulicProperties(int waterRetentionCurve,
  void DLL_EXPORT __STDCALL initializeBalance()
 {
     InitializeBalanceWater();
+
     if (myStructure.computeHeat)
+    {
         initializeBalanceHeat();
+    }
     else
+    {
         balanceWholePeriod.heatMBR = 1.;
+    }
 }
 
 

agrolib/soilFluxes3D/water.cpp

@@ -323,7 +323,7 @@ bool waterFlowComputation(double deltaT)
         isValidStep = waterBalance(deltaT, approximationNr);
         if (getForcedHalvedTime()) return (false);
         }
-    while ((!isValidStep) && (++approximationNr < myParameters.maxApproximationsNumber));
+    while ((! isValidStep) && (++approximationNr < myParameters.maxApproximationsNumber));
 
     return isValidStep;
  }
@@ -368,7 +368,7 @@ bool computeWater(double maxTime, double *acceptedTime)
 
         isStepOK = waterFlowComputation(*acceptedTime);
 
-        if (!isStepOK) restoreWater();
+        if (! isStepOK) restoreWater();
     }
     return (isStepOK);
 }