Skip to content

Commit

Permalink
Merge commit 'dd99bea0a799b4df858991f7fe576e146679da8d'
Browse files Browse the repository at this point in the history
  • Loading branch information
ftomei committed Jun 17, 2024
2 parents a4fb329 + dd99bea commit 0e41fa8
Show file tree
Hide file tree
Showing 3 changed files with 85 additions and 14 deletions.
15 changes: 8 additions & 7 deletions agrolib/interpolation/interpolation.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -933,7 +933,7 @@ float modifiedShepardIdw(vector <Crit3DInterpolationDataPoint> &myPoints,
{
t[i] = 0;
for (j=0; j < validPoints.size(); j++)
if (i != j && S[i] > 0)
if (i != j && S[i] > 0 && validPoints[i].distance > 0 && validPoints[j].distance > 0)
{
cosine = ((X - (float)validPoints[i].point->utm.x) * (X - (float)validPoints[j].point->utm.x) + (Y - (float)validPoints[i].point->utm.y) * (Y - (float)validPoints[j].point->utm.y)) / (validPoints[i].distance * validPoints[j].distance);
t[i] = t[i] + S[j] * (1 - cosine);
Expand Down Expand Up @@ -1060,8 +1060,9 @@ void localSelection(vector <Crit3DInterpolationDataPoint> &inputPoints, vector <
r1 += stepRadius;
}

for (i=0; i< selectedPoints.size(); i++)
selectedPoints[i].regressionWeight = MAXVALUE(1 - selectedPoints[i].distance / maxDistance,2*EPSILON);
if (maxDistance != 0)
for (i=0; i< selectedPoints.size(); i++)
selectedPoints[i].regressionWeight = MAXVALUE(1 - selectedPoints[i].distance / (maxDistance*maxDistance*maxDistance),EPSILON);

mySettings.setLocalRadius(maxDistance);
}
Expand Down Expand Up @@ -1357,7 +1358,7 @@ bool setAllFittingRanges(Crit3DProxyCombination myCombination, Crit3DInterpolati
{
mySettings->getProxy(i)->setFittingFunctionName(piecewiseTwo);
if ((mySettings->getProxy(i)->getFittingParametersRange().empty()))
tempParam = {-200, min-2, 0, -0.006, 1800, max+2, 0.01, 0};
tempParam = {-200, min-2, 0.001, -0.006, 1800, max+2, 0.01, 0.0015};
else
{
tempParam = mySettings->getProxy(i)->getFittingParametersRange();
Expand All @@ -1369,7 +1370,7 @@ bool setAllFittingRanges(Crit3DProxyCombination myCombination, Crit3DInterpolati
{
mySettings->getProxy(i)->setFittingFunctionName(piecewiseThreeFree);
if ((mySettings->getProxy(i)->getFittingParametersRange().empty()))
tempParam = {-200, min-2, 100, 0.001, -0.006, -0.006, 1800, max+2, 1000, 0.007, 0, 0};
tempParam = {-200, min-2, 100, 0.001, -0.006, -0.006, 1800, max+2, 1000, 0.007, 0.0015, 0.0015};
else
{
tempParam = mySettings->getProxy(i)->getFittingParametersRange();
Expand All @@ -1381,7 +1382,7 @@ bool setAllFittingRanges(Crit3DProxyCombination myCombination, Crit3DInterpolati
{
mySettings->getProxy(i)->setFittingFunctionName(piecewiseThree);
if ((mySettings->getProxy(i)->getFittingParametersRange().empty()))
tempParam = {-200, min-2, 100, 0.001, -0.006, 1800, max+2, 1000, 0.007, 0};
tempParam = {-200, min-2, 100, 0.001, -0.006, 1800, max+2, 1000, 0.007, 0.0015};
else
{
tempParam = mySettings->getProxy(i)->getFittingParametersRange();
Expand Down Expand Up @@ -1803,7 +1804,7 @@ bool multipleDetrendingElevation(Crit3DProxyCombination elevationCombination, st
//std::vector <std::vector<double>> parameters;
std::vector<std::function<double(double, std::vector<double>&)>> myFunc;

unsigned int nrMaxStep = 20;
unsigned int nrMaxStep = 100;
if (parameters.empty())
nrMaxStep *= 10;

Expand Down
82 changes: 77 additions & 5 deletions agrolib/mathFunctions/furtherMathFunctions.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1971,7 +1971,7 @@ namespace interpolation
std::vector<double> ySim(nrData);

int counter = 0;
//srand (unsigned(time(nullptr)));
/*srand (unsigned(time(nullptr)));
std::random_device rd;
std::mt19937 gen(rd());
std::normal_distribution<double> normal_dis(0.5, 0.5);
Expand All @@ -1989,9 +1989,80 @@ namespace interpolation
discreteParameters.push_back(tempDiscrete);
tempDiscrete.clear();
}*/

//grigliato

std::vector<double> steps;
if (parameters.size() == 4)
steps = {50.0, 0.5, 0.00075, 0.00075};
else if (parameters.size() == 6)
steps = {50.0, 0.5, 50.0, 0.00075, 0.00075, 0.00075};
else return false;

const int numSteps = 40;
int directions[] = {1, -1};
int numParamsToVary = parameters.size();
std::vector<double> firstGuessParam = parameters;

for (int step = 1; step <= numSteps; ++step)
{
for (int dir = 0; dir < 2; ++dir)
{
for (int paramIndex = 0; paramIndex < numParamsToVary; ++paramIndex)
{

fittingMarquardt_nDimension_noSquares_singleFunction(func,parametersMin,
parametersMax,parameters,
parametersDelta,maxIterationsNr,
myEpsilon,x,y,weights);

for (i=0;i<nrData;i++)
{
ySim[i]= func(x[i], parameters);
}
R2 = computeWeighted_R2(y,ySim,weights);

if (R2 > (bestR2-deltaR2))
{
for (j=0;j<nrMinima-1;j++)
{
R2Previous[j] = R2Previous[j+1];
}
if (R2 > (bestR2))
{
for (j=0; j<nrParameters; j++)
{
bestParameters[j] = parameters[j];
}
bestR2 = R2;
}
R2Previous[nrMinima-1] = R2;

for (j=0; j<nrParameters; j++)
{
bestParameters[j] = parameters[j];
}
}
counter++;

if (dir == 0)
parameters[paramIndex] = MINVALUE(firstGuessParam[paramIndex] + directions[dir] * step * steps[paramIndex], parametersMax[paramIndex]);
else
parameters[paramIndex] = MAXVALUE(firstGuessParam[paramIndex] + directions[dir] * step * steps[paramIndex], parametersMin[paramIndex]);

}
}

if ((counter > nrTrials) || ((R2Previous[0] != NODATA) && fabs(R2Previous[0]-R2Previous[nrMinima-1]) < deltaR2 ))
break;


}
//end grigliato

do
//random
/*do
{
fittingMarquardt_nDimension_noSquares_singleFunction(func,parametersMin,
parametersMax,parameters,
Expand Down Expand Up @@ -2030,13 +2101,13 @@ namespace interpolation
}
counter++;
/*for (j=0; j<nrParameters; j++)
for (j=0; j<nrParameters; j++)
{
do {
truncNormal = normal_dis(gen);
} while(truncNormal <= 0.0 || truncNormal >= 1.0);
parameters[j] = parametersMin[j] + (truncNormal)*(parametersMax[j]-parametersMin[j]);
}*/
}
int indice = 0;
for (j=0; j<nrParameters-2;j++)
Expand All @@ -2046,8 +2117,9 @@ namespace interpolation
} while (indice < 0 || indice >= discreteParameters[j].size());
parameters[j] = discreteParameters[j][indice];
}
} while( (counter < nrTrials) && !(R2Previous[0] > 0.797 && R2Previous[nrMinima-1] > 0.8) && ((R2Previous[0] == NODATA) || fabs(R2Previous[0]-R2Previous[nrMinima-1]) > deltaR2 ));
} while( (counter < nrTrials) && !(R2Previous[0] > 0.797 && R2Previous[nrMinima-1] > 0.8) && ((R2Previous[0] == NODATA) || fabs(R2Previous[0]-R2Previous[nrMinima-1]) > deltaR2 ));
*/
for (j=0; j<nrParameters; j++)
{
parameters[j] = bestParameters[j];
Expand Down
2 changes: 0 additions & 2 deletions agrolib/project/project.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -2443,10 +2443,8 @@ bool Project::interpolationDemLocalDetrending(meteoVariable myVar, const Crit3DT
return false;
}


myRaster->value[row][col] = interpolate(subsetInterpolationPoints, &interpolationSettings, meteoSettings,
myVar, x, y, z, proxyValues, true);

myRaster->setParametersForRowCol(row, col, interpolationSettings.getFittingParameters());
interpolationSettings.setCurrentCombination(interpolationSettings.getSelectedCombination());
}
Expand Down

0 comments on commit 0e41fa8

Please sign in to comment.