Skip to content

Commit

Permalink
update vine3D
Browse files Browse the repository at this point in the history
  • Loading branch information
@ftomei
ftomei committed Jul 26, 2024
1 parent a31fd88 commit 4687a03
Show file tree
Hide file tree
Showing 16 changed files with 213 additions and 191 deletions.
2 changes: 1 addition & 1 deletion DATA/PROJECT/VINE3D_test/SETTINGS/parameters.ini
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ relhum_tolerance=102
water_table_maximum_depth=300

[meteo]
min_percentage=60
min_percentage=80
prec_threshold=0.2
samani_coefficient=0.17
thom_threshold=24
Expand Down
10 changes: 0 additions & 10 deletions agrolib/crop/landUnit.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -59,13 +59,3 @@ bool loadLandUnitList(const QSqlDatabase &dbCrop, std::vector<Crit3DLandUnit> &l

return true;
}


int getLandUnitIndex(const std::vector<Crit3DLandUnit> &landUnitList, int idLandUnit)
{
for (int index = 0; index < int(landUnitList.size()); index++)
if (landUnitList[index].id == idLandUnit)
return index;

return NODATA;
}
2 changes: 0 additions & 2 deletions agrolib/crop/landUnit.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,5 @@

bool loadLandUnitList(const QSqlDatabase &dbCrop, std::vector<Crit3DLandUnit> &landUnitList, QString &errorStr);

int getLandUnitIndex(const std::vector<Crit3DLandUnit> &landUnitList, int idLandUnit);


#endif // LANDUNIT_H
115 changes: 61 additions & 54 deletions agrolib/grapevine/grapevine.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,9 +21,7 @@ const float LAIMIN = 0.01f;


Vine3D_Grapevine::Vine3D_Grapevine()
{
}

{}

bool Vine3D_Grapevine::compute(bool computeDaily, int secondsPerStep, Crit3DModelCase* modelCase, double chlorophyll)
{
Expand Down Expand Up @@ -106,9 +104,9 @@ bool Vine3D_Grapevine::compute(bool computeDaily, int secondsPerStep, Crit3DMode
grassTranspiration(modelCase);

return(true);

}


void Vine3D_Grapevine::resetLayers()
{
for (int i=0 ; i < nrMaxLayers ; i++)
Expand All @@ -121,6 +119,7 @@ void Vine3D_Grapevine::resetLayers()
}
}


void Vine3D_Grapevine::initializeLayers(int myMaxLayers)
{
nrMaxLayers = myMaxLayers;
Expand Down Expand Up @@ -156,22 +155,25 @@ bool Vine3D_Grapevine::setWeather(double meanDailyTemp, double temp, double irra
myMeanDailyTemperature = meanDailyTemp;
double deltaRelHum = MAXVALUE(100.0 - myRelativeHumidity, 0.01);
myVaporPressureDeficit = 0.01 * deltaRelHum * 613.75 * exp(17.502 * myInstantTemp / (240.97 + myInstantTemp));
//globalRadiation = globRad;
if ((int(prec) != NODATA) && (int(temp) != NODATA) && (int(windSpeed) != NODATA) && (int(irradiance) != NODATA) && (int(relativeHumidity) != NODATA) && (int(atmosphericPressure) != NODATA)) isReadingOK = true ;
return isReadingOK ;

if ((int(prec) != NODATA) && (int(temp) != NODATA) && (int(windSpeed) != NODATA)
&& (int(irradiance) != NODATA) && (int(relativeHumidity) != NODATA) && (int(atmosphericPressure) != NODATA))
isReadingOK = true;

return isReadingOK;
}


bool Vine3D_Grapevine::setDerivedVariables(double diffuseIrradiance, double directIrradiance, double cloudIndex, double sunElevation)
{
bool isReadingOK = false;
myDiffuseIrradiance = diffuseIrradiance ;
myDirectIrradiance = directIrradiance ;
//myLongWaveIrradiance = longWaveIrradiance ;
myCloudiness = MINVALUE(1,MAXVALUE(0,cloudIndex)) ;
//myAirVapourPressure = airVapourPressure ;
mySunElevation = sunElevation;
if (int(sunElevation) != NODATA && int(diffuseIrradiance) != NODATA && int(directIrradiance) != NODATA
&& int(cloudIndex) != NODATA) isReadingOK = true ;
if (int(sunElevation) != NODATA && int(diffuseIrradiance) != NODATA && int(directIrradiance) != NODATA && int(cloudIndex) != NODATA)
isReadingOK = true ;

return isReadingOK ;
}

Expand All @@ -185,32 +187,45 @@ void Vine3D_Grapevine::initializeWaterStress(Crit3DModelCase* modelCase)
}


bool Vine3D_Grapevine::setSoilProfile(Crit3DModelCase* modelCase, double* myWiltingPoint, double* myFieldCapacity,
double* myPsiSoilProfile, double* mySoilWaterContentProfile,
double* mySoilWaterContentFC, double* mySoilWaterContentWP)
bool Vine3D_Grapevine::setSoilProfile(Crit3DModelCase* modelCase, std::vector<double>& myWiltingPoint, std::vector<double>& myFieldCapacity,
std::vector<double>& myPsiSoilProfile, std::vector<double>& mySoilWaterContentProfile,
std::vector<double>& mySoilWaterContentFC, std::vector<double>& mySoilWaterContentWP)
{
double psiSoilProfile;
double logPsiSoilAverage = 0.;
double logPsiFCAverage = 0.;
double soilFieldCapacity;
// check last Soil layer
int lastLayer = modelCase->soilLayersNr - 1;
for (int i = 1; i < modelCase->soilLayersNr; i++)
{
if (isEqual(myPsiSoilProfile[i], NODATA) || isEqual(myFieldCapacity[i], NODATA) || isEqual(mySoilWaterContentProfile[i], NODATA)
|| isEqual(mySoilWaterContentProfile[i], NODATA) || isEqual(mySoilWaterContentWP[i], NODATA))
{
lastLayer = i-1;
break;
}
}
if (lastLayer == 0)
return false;

psiSoilAverage = 0.;
psiFieldCapacityAverage = 0.;

if (int(myWiltingPoint[int(modelCase->soilLayersNr / 2)]) == NODATA)
return false;
// initialize vector
for (int i = 1; i < modelCase->soilLayersNr; i++)
{
transpirationLayer[i] = 0.;
transpirationCumulatedGrass[i] = 0.;
fractionTranspirableSoilWaterProfile[i] = 0.;
}

wiltingPoint = myWiltingPoint[int(modelCase->soilLayersNr / 2)] / 101.97; // conversion from mH2O to MPa
int midLayer = (lastLayer + 1) / 2;
wiltingPoint = myWiltingPoint[midLayer] / 101.97; // conversion from mH2O to MPa

//layer 0: surface, no soil
for (int i = 1; i < modelCase->soilLayersNr; i++)
double logPsiSoilAverage = 0.;
double logPsiFCAverage = 0.;
for (int i = 1; i <= lastLayer; i++)
{
if (isEqual(myPsiSoilProfile[i], NODATA) || isEqual(myFieldCapacity[i], NODATA) || isEqual(mySoilWaterContentProfile[i], NODATA)
|| isEqual(mySoilWaterContentProfile[i], NODATA) || isEqual(mySoilWaterContentWP[i], NODATA))
return false;

soilFieldCapacity = myFieldCapacity[i]/101.97; // conversion from mH2O to MPa
psiSoilProfile = MINVALUE(myPsiSoilProfile[i],-1.)/101.97 ; // conversion from mH2O to MPa
double soilFieldCapacity = myFieldCapacity[i] / 101.97; // conversion from mH2O to MPa
double psiSoilProfile = MINVALUE(myPsiSoilProfile[i], -1.) / 101.97 ; // conversion from mH2O to MPa
logPsiSoilAverage += log(-psiSoilProfile) * modelCase->rootDensity[i];
logPsiFCAverage += log(-soilFieldCapacity) * modelCase->rootDensity[i];
}
Expand All @@ -219,46 +234,35 @@ bool Vine3D_Grapevine::setSoilProfile(Crit3DModelCase* modelCase, double* myWilt
psiFieldCapacityAverage = -exp(logPsiFCAverage);
fractionTranspirableSoilWaterAverage = 0;

double waterContent, waterContentFC, waterContentWP;

for (int i = 0; i < modelCase->soilLayersNr; i++)
for (int i = 0; i <= lastLayer; i++)
{
waterContent = mySoilWaterContentProfile[i];
waterContentFC = mySoilWaterContentFC[i];
waterContentWP = mySoilWaterContentWP[i];
double waterContent = mySoilWaterContentProfile[i];
double waterContentFC = mySoilWaterContentFC[i];
double waterContentWP = mySoilWaterContentWP[i];

fractionTranspirableSoilWaterProfile[i] = MAXVALUE(0, MINVALUE(1, (waterContent - waterContentWP) / (waterContentFC - waterContentWP)));
fractionTranspirableSoilWaterAverage += fractionTranspirableSoilWaterProfile[i] * modelCase->rootDensity[i];
transpirationLayer[i] = 0.;
transpirationCumulatedGrass[i] = 0. ;
}

return true ;
}

bool Vine3D_Grapevine::setStatePlant(TstatePlant myStatePlant, bool isVineyard)

bool Vine3D_Grapevine::setStatePlant(TstatePlant myStatePlant, bool isVineyard_)
{
statePlant = myStatePlant;
this->statePlant.outputPlant.transpirationNoStress = 0.;

if (! isVineyard)
if (! isVineyard_)
{
statePlant.outputPlant.brixBerry = NODATA;
statePlant.statePheno.stage = NODATA;
statePlant.stateGrowth.cumulatedBiomass = NODATA;
statePlant.stateGrowth.leafAreaIndex = NODATA;
statePlant.stateGrowth.fruitBiomass = NODATA;
}
return true;
}

TstatePlant Vine3D_Grapevine::getStatePlant()
{
return(this->statePlant);
}

ToutputPlant Vine3D_Grapevine::getOutputPlant()
{
return (this->statePlant.outputPlant);
return true;
}


Expand All @@ -268,13 +272,13 @@ void Vine3D_Grapevine::getFixSimulationParameters()
parameterBindiMigliettaFix.a = -0.28;
parameterBindiMigliettaFix.b = 0.04;
parameterBindiMigliettaFix.c = -0.015;
//parameterBindiMigliettaFix.baseTemperature = 10;
//parameterBindiMigliettaFix.tempMaxThreshold = 35;
parameterBindiMigliettaFix.extinctionCoefficient = 0.5;
parameterBindiMigliettaFix.shadedSurface = 0.8;

// Wang Leuning parameters
parameterWangLeuningFix.stomatalConductanceMin = 0.008;
parameterWangLeuningFix.optimalTemperatureForPhotosynthesis = 298.15;

// fenovitis parameters
parameterPhenoVitisFix.a= 0.005;
parameterPhenoVitisFix.optimalChillingTemp = 2.8;
Expand Down Expand Up @@ -525,6 +529,8 @@ void Vine3D_Grapevine::radiationAbsorption()
{
sunlit.leafAreaIndex = 0.0 ;
sunlit.absorbedPAR = 0.0 ;

// TODO: non servono?
sunlitAbsorbedNIR = 0.0 ;
sunlitAbsorbedLW = 0.0 ;
sunlit.isothermalNetRadiation = 0.0 ;
Expand All @@ -537,12 +543,13 @@ void Vine3D_Grapevine::radiationAbsorption()
shadedAbsorbedLW= dum[16] * (UPSCALINGFUNC(diffuseLightK,statePlant.stateGrowth.leafAreaIndex) - UPSCALINGFUNC(directLightK+diffuseLightK,statePlant.stateGrowth.leafAreaIndex)) ;
shaded.isothermalNetRadiation = shaded.absorbedPAR + shadedAbsorbedNIR + shadedAbsorbedLW ;
}

// Convert absorbed PAR into units of mol m-2 s-1
sunlit.absorbedPAR *= 4.57E-6 ;
shaded.absorbedPAR *= 4.57E-6 ;

}


void Vine3D_Grapevine::leafTemperature()
{
// taken from Hydrall Model, Magnani UNIBO
Expand Down Expand Up @@ -695,7 +702,8 @@ void Vine3D_Grapevine::aerodynamicalCoupling()
{
sunlitDeltaTemp = 0.0;
}
sunlitDeltaTemp = 0.0;
sunlitDeltaTemp = 0.0; // TODO: check

sunlit.leafTemperature = myInstantTemp + sunlitDeltaTemp + ZEROCELSIUS ; //sunlit big-leaf
stomatalConductanceWater = 10.0/shaded.leafAreaIndex ; //dummy stom res for shaded big-leaf
//if (shaded.isothermalNetRadiation > 100) stomatalConductanceWater *= pow(100/shaded.isothermalNetRadiation,0.5);
Expand Down Expand Up @@ -802,10 +810,9 @@ void Vine3D_Grapevine::upscale(TVineCultivar *cultivar)
sunlit.darkRespiration = 0.0;
sunlit.maximalCarboxylationRate = 0.0;
}


}


void Vine3D_Grapevine::photosynthesisKernel(TVineCultivar* cultivar, double COMP,double GAC,double GHR,double GSCD,double J,double KC,double KO
,double RD,double RNI,double STOMWL,double VCmax,double *ASS,double *GSC,double *TR)
{
Expand Down
29 changes: 14 additions & 15 deletions agrolib/grapevine/grapevine.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,9 +41,9 @@
enum phenoStage {endoDormancy, ecoDormancy , budBurst , flowering , fruitSet, veraison, physiologicalMaturity, vineSenescence};
enum TfieldOperation {irrigationOperation, grassSowing, grassRemoving, trimming, leafRemoval,
clusterThinning, harvesting, tartaricAnalysis};
enum Crit3DLanduse {landuse_nodata, landuse_bare, landuse_vineyard};
enum GrapevineLanduse {landuse_nodata, landuse_bare, landuse_vineyard};

const std::map<std::string, Crit3DLanduse> landuseNames = {
const std::map<std::string, GrapevineLanduse> landuseNames = {
{ "UNDEFINED", landuse_nodata },
{ "BARESOIL", landuse_bare },
{ "VINEYARD", landuse_vineyard}
Expand Down Expand Up @@ -191,7 +191,7 @@

struct Crit3DModelCase {
int id;
Crit3DLanduse landuse;
GrapevineLanduse landuse;
int soilIndex;

float shootsPerPlant;
Expand Down Expand Up @@ -300,8 +300,7 @@
double psiSoilAverage;
double psiFieldCapacityAverage;

//double* layerRootDensity;
double totalStomatalConductance, totalStomatalConductanceNoStress ;
double totalStomatalConductance, totalStomatalConductanceNoStress;
double transpirationInstant;
double* currentProfile;
double* transpirationInstantLayer; //molH2O m^-2 s^-1
Expand All @@ -327,8 +326,6 @@
bool isAmphystomatic ;
double specificLeafArea ;
double alphaLeuning ;
//double leafNitrogen ;
//double entropicFactorCarboxyliation,entropicFactorElectronTransporRate ;

Vine3D_SunShade shaded ;
Vine3D_SunShade sunlit ;
Expand Down Expand Up @@ -374,7 +371,7 @@

double getWaterStressByPsiSoil(double myPsiSoil,double psiSoilStressParameter,double exponentialFactorForPsiRatio);
double getWaterStressSawFunction(int index, TVineCultivar *cultivar);
//bool getExtractedWaterFromGrassTranspirationandEvaporation(double* myWaterExtractionProfile);

double getWaterStressSawFunctionAverage(TVineCultivar* cultivar);
double getGrassTranspiration(double stress, double laiGrassMax, double sensitivityToVPD, double fieldCoverByPlant);
double getFallowTranspiration(double stress, double laiGrassMax, double sensitivityToVPD);
Expand All @@ -392,6 +389,11 @@
public:
Vine3D_Grapevine();

TstatePlant getStatePlant(){ return statePlant; }
ToutputPlant getOutputPlant() { return statePlant.outputPlant; }

bool setStatePlant(TstatePlant myStatePlant, bool isVineyard_);

//void initializeGrapevineModel(TVineCultivar* cultivar, double secondsPerStep);
void initializeLayers(int myMaxLayers);
bool initializeStatePlant(int doy, Crit3DModelCase *vineField);
Expand All @@ -409,15 +411,12 @@
double prec , double relativeHumidity , double windSpeed, double atmosphericPressure);
bool setDerivedVariables (double diffuseIrradiance, double directIrradiance,
double cloudIndex, double sunElevation);
bool setSoilProfile(Crit3DModelCase *modelCase, double* myWiltingPoint, double *myFieldCapacity,
double *myPsiSoilProfile , double *mySoilWaterContentProfile,
double* mySoilWaterContentFC, double* mySoilWaterContentWP);
bool setStatePlant(TstatePlant myStatePlant, bool isVineyard);

TstatePlant getStatePlant();
ToutputPlant getOutputPlant();
bool setSoilProfile(Crit3DModelCase* modelCase, std::vector<double>& myWiltingPoint, std::vector<double>& myFieldCapacity,
std::vector<double>& myPsiSoilProfile, std::vector<double>& mySoilWaterContentProfile,
std::vector<double>& mySoilWaterContentFC, std::vector<double>& mySoilWaterContentWP);

double* getExtractedWater(Crit3DModelCase* modelCase);
//bool getOutputPlant(int hour, ToutputPlant *outputPlant);
double getStressCoefficient();
double getRealTranspirationGrapevine(Crit3DModelCase *modelCase);
double getRealTranspirationGrass(Crit3DModelCase *modelCase);
Expand Down
9 changes: 7 additions & 2 deletions agrolib/project/project.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -934,8 +934,13 @@ Crit3DTime Project::getCrit3DCurrentTime()
QDateTime Project::getCurrentTime()
{
QDateTime myDateTime;
myDateTime.setDate(this->currentDate);
return myDateTime.addSecs(this->currentHour * HOUR_SECONDS);
if (gisSettings.isUTC)
{
myDateTime.setTimeSpec(Qt::UTC);
}

myDateTime.setDate(currentDate);
return myDateTime.addSecs(currentHour * HOUR_SECONDS);
}


Expand Down
Loading

0 comments on commit 4687a03

Please sign in to comment.