Skip to content

Commit

Permalink
Merge branch 'master' of https://github.com/ARPA-SIMC/CRITERIA1D
Browse files Browse the repository at this point in the history
  • Loading branch information
@gantolini
gantolini committed Jan 29, 2024
2 parents 6cb80fe + 7121591 commit cd01594
Show file tree
Hide file tree
Showing 33 changed files with 1,041 additions and 178 deletions.
Binary file modified DOC/CRITERIA1D_technical_manual.pdf
View file
Binary file not shown.
40 changes: 34 additions & 6 deletions agrolib/criteriaModel/criteria1DCase.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -361,7 +361,7 @@ bool Crit1DCase::computeNumericalFluxes(const Crit3DDate &myDate, std::string &e
}


void Crit1DCase::saveWaterContent()
void Crit1DCase::storeWaterContent()
{
prevWaterContent.clear();
prevWaterContent.resize(soilLayers.size());
Expand Down Expand Up @@ -537,7 +537,7 @@ bool Crit1DCase::computeDailyModel(Crit3DDate &myDate, std::string &error)

output.dailyIrrigation = 0;
// Water fluxes (first computation)
saveWaterContent();
storeWaterContent();
if (! computeWaterFluxes(myDate, error)) return false;
// Irrigation
double irrigation = checkIrrigationDemand(doy, prec, precTomorrow, output.dailyMaxTranspiration);
Expand Down Expand Up @@ -653,11 +653,11 @@ double Crit1DCase::getTotalWaterContent()


/*!
* \brief getWaterContent
* \brief getVolumetricWaterContent
* \param computationDepth = computation soil depth [cm]
* \return volumetric water content at specific depth [-]
* \return volumetric water content (soil moisture) at specific depth [m3 m-3]
*/
double Crit1DCase::getWaterContent(double computationDepth)
double Crit1DCase::getVolumetricWaterContent(double computationDepth)
{
computationDepth /= 100; // [cm] --> [m]

Expand All @@ -679,6 +679,34 @@ double Crit1DCase::getWaterContent(double computationDepth)
}


/*!
* \brief getDegreeOfSaturation
* \param computationDepth = computation soil depth [cm]
* \return degree of saturation at specific depth [-]
*/
double Crit1DCase::getDegreeOfSaturation(double computationDepth)
{
computationDepth /= 100; // [cm] --> [m]
if (computationDepth <= 0 || computationDepth > mySoil.totalDepth)
{
return NODATA;
}

double upperDepth, lowerDepth;
for (unsigned int i = 1; i < soilLayers.size(); i++)
{
upperDepth = soilLayers[i].depth - soilLayers[i].thickness * 0.5;
lowerDepth = soilLayers[i].depth + soilLayers[i].thickness * 0.5;
if (computationDepth >= upperDepth && computationDepth <= lowerDepth)
{
return soilLayers[i].waterContent / soilLayers[i].SAT;
}
}

return NODATA;
}


/*!
* \brief getWaterPotential
* \param computationDepth = computation soil depth [cm]
Expand Down Expand Up @@ -778,7 +806,7 @@ double Crit1DCase::getWaterDeficitSum(double computationDepth)

/*!
* \brief getWaterCapacity
* \param computationDepth = computation soil depth [cm]
* \param computationDepth = computation soil depth [cm]
* \return sum of available water capacity (FC-WP) from zero to computationDepth (mm)
*/
double Crit1DCase::getWaterCapacitySum(double computationDepth)
Expand Down
5 changes: 3 additions & 2 deletions agrolib/criteriaModel/criteria1DCase.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,8 @@
bool initializeWaterContent(Crit3DDate myDate);
bool computeDailyModel(Crit3DDate &myDate, std::string &error);

double getWaterContent(double computationDepth);
double getVolumetricWaterContent(double computationDepth);
double getDegreeOfSaturation(double computationDepth);
double getWaterPotential(double computationDepth);
double getFractionAW(double computationDepth);
double getSlopeStability(double computationDepth);
Expand All @@ -99,7 +100,7 @@
bool computeNumericalFluxes(const Crit3DDate &myDate, std::string &error);
bool computeWaterFluxes(const Crit3DDate &myDate, std::string &error);
double checkIrrigationDemand(int doy, double currentPrec, double nextPrec, double maxTranspiration);
void saveWaterContent();
void storeWaterContent();
void restoreWaterContent();
double getTotalWaterContent();

Expand Down
28 changes: 25 additions & 3 deletions agrolib/criteriaModel/criteria1DProject.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,9 +63,11 @@ void Crit1DProject::initialize()
lastSimulationDate = QDate(1800,1,1);

outputString = "";

// specific outputs
waterDeficitDepth.clear();
waterContentDepth.clear();
degreeOfSaturationDepth.clear();
waterPotentialDepth.clear();
availableWaterDepth.clear();
fractionAvailableWaterDepth.clear();
Expand Down Expand Up @@ -257,6 +259,12 @@ bool Crit1DProject::readSettings()
projectError = "Wrong water content depth in " + configFileName;
return false;
}
depthList = projectSettings->value("degreeOfSaturation").toStringList();
if (! setVariableDepth(depthList, degreeOfSaturationDepth))
{
projectError = "Wrong degree of saturation depth in " + configFileName;
return false;
}
depthList = projectSettings->value("waterPotential").toStringList();
if (! setVariableDepth(depthList, waterPotentialDepth))
{
Expand Down Expand Up @@ -1544,7 +1552,12 @@ bool Crit1DProject::createOutputTable(QString &myError)
// specific depth variables
for (unsigned int i = 0; i < waterContentDepth.size(); i++)
{
QString fieldName = "SWC_" + QString::number(waterContentDepth[i]);
QString fieldName = "VWC_" + QString::number(waterContentDepth[i]);
queryString += ", " + fieldName + " REAL";
}
for (unsigned int i = 0; i < degreeOfSaturationDepth.size(); i++)
{
QString fieldName = "DEGSAT_" + QString::number(degreeOfSaturationDepth[i]);
queryString += ", " + fieldName + " REAL";
}
for (unsigned int i = 0; i < waterPotentialDepth.size(); i++)
Expand Down Expand Up @@ -1604,7 +1617,12 @@ void Crit1DProject::updateOutput(Crit3DDate myDate, bool isFirst)
// specific depth variables
for (unsigned int i = 0; i < waterContentDepth.size(); i++)
{
QString fieldName = "SWC_" + QString::number(waterContentDepth[i]);
QString fieldName = "VWC_" + QString::number(waterContentDepth[i]);
outputString += ", " + fieldName;
}
for (unsigned int i = 0; i < degreeOfSaturationDepth.size(); i++)
{
QString fieldName = "DEGSAT_" + QString::number(degreeOfSaturationDepth[i]);
outputString += ", " + fieldName;
}
for (unsigned int i = 0; i < waterPotentialDepth.size(); i++)
Expand Down Expand Up @@ -1669,7 +1687,11 @@ void Crit1DProject::updateOutput(Crit3DDate myDate, bool isFirst)
// specific depth variables
for (unsigned int i = 0; i < waterContentDepth.size(); i++)
{
outputString += "," + QString::number(myCase.getWaterContent(waterContentDepth[i]), 'g', 4);
outputString += "," + QString::number(myCase.getVolumetricWaterContent(waterContentDepth[i]), 'g', 4);
}
for (unsigned int i = 0; i < degreeOfSaturationDepth.size(); i++)
{
outputString += "," + QString::number(myCase.getDegreeOfSaturation(degreeOfSaturationDepth[i]), 'g', 4);
}
for (unsigned int i = 0; i < waterPotentialDepth.size(); i++)
{
Expand Down
1 change: 1 addition & 0 deletions agrolib/criteriaModel/criteria1DProject.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -91,6 +91,7 @@

// specific output
std::vector<int> waterContentDepth;
std::vector<int> degreeOfSaturationDepth;
std::vector<int> waterPotentialDepth;
std::vector<int> waterDeficitDepth;
std::vector<int> awcDepth;
Expand Down
2 changes: 2 additions & 0 deletions agrolib/crop/root.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -591,6 +591,7 @@ namespace root
}
}

// normalize root density
if (rootDensitySumSubset != rootDensitySum)
{
double ratio = rootDensitySum / rootDensitySumSubset;
Expand All @@ -600,6 +601,7 @@ namespace root
}
}

// find first and last root layers
myCrop->roots.firstRootLayer = 0;
unsigned int layer = 0;
while (layer < nrLayers && myCrop->roots.rootDensity[layer] == 0.0)
Expand Down
118 changes: 47 additions & 71 deletions agrolib/dbMeteoGrid/dbMeteoGrid.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
#include "meteoGrid.h"
#include "basicMath.h"
#include "utilities.h"
#include "meteoPoint.h"
#include "commonConstants.h"

#include <iostream>
Expand Down Expand Up @@ -54,7 +55,6 @@ bool Crit3DMeteoGridDbHandler::parseXMLFile(QString xmlFileName, QDomDocument* x

bool Crit3DMeteoGridDbHandler::parseXMLGrid(QString xmlFileName, QString *myError)
{

QDomDocument xmlDoc;

if (! parseXMLFile(xmlFileName, &xmlDoc, myError)) return false;
Expand Down Expand Up @@ -1702,15 +1702,15 @@ bool Crit3DMeteoGridDbHandler::updateMeteoGridDate(QString &myError)
}


bool Crit3DMeteoGridDbHandler::loadGridDailyData(QString &myError, const QString &meteoPoint, const QDate &firstDate, const QDate &lastDate)
bool Crit3DMeteoGridDbHandler::loadGridDailyData(QString &myError, const QString &meteoPointId, const QDate &firstDate, const QDate &lastDate)
{
myError = "";
QString tableD = _tableDaily.prefix + meteoPoint + _tableDaily.postFix;
QString tableD = _tableDaily.prefix + meteoPointId + _tableDaily.postFix;

unsigned row, col;
if ( !_meteoGrid->findMeteoPointFromId(&row, &col, meteoPoint.toStdString()) )
if ( !_meteoGrid->findMeteoPointFromId(&row, &col, meteoPointId.toStdString()) )
{
myError = "Missing MeteoPoint id";
myError = "Missing meteoPoint id: " + meteoPointId;
return false;
}

Expand All @@ -1723,6 +1723,7 @@ bool Crit3DMeteoGridDbHandler::loadGridDailyData(QString &myError, const QString
{
if (firstDate > _lastDailyDate || lastDate < _firstDailyDate)
{
myError = "Missing data in this time interval.";
return false;
}
}
Expand Down Expand Up @@ -3631,36 +3632,42 @@ bool Crit3DMeteoGridDbHandler::saveLogProcedures(QString *myError, QString nameP
/*!
* \brief ExportDailyDataCsv
* export gridded daily meteo data to csv files
* \param isTPrec save only variables: Tmin, Tmax, Tavg, Prec
* \param idListFileName filename of cells id list (list by columns)
* if idListFile == "" save ALL cells
* \param variableList list of meteo variables
* \param idListFileName text file of cells id list by columns - if idListFileName is empty save ALL cells
* \param outputPath path for output files
* \return true on success, false otherwise
*/
bool Crit3DMeteoGridDbHandler::exportDailyDataCsv(QString &errorStr, bool isTPrec, QDate firstDate, QDate lastDate,
QString idListFileName, QString outputPath)
bool Crit3DMeteoGridDbHandler::exportDailyDataCsv(QString &errorStr, QList<meteoVariable> variableList,
QDate firstDate, QDate lastDate, QString idListFileName, QString outputPath)
{
errorStr = "";

// check output path
if (outputPath == "")
{
errorStr = "Missing output path.";
return false;
}

QDir outDir(outputPath);
// make directory
if (! outDir.exists())
{
if (! outDir.mkpath(outputPath))
{
errorStr = "Wrong outputPath, unable to create this directory: " + outputPath;
errorStr = "Wrong output path, unable to create directory: " + outputPath;
return false;
}
}
outputPath = outDir.absolutePath();

bool isList = (idListFileName != "");
bool isList = (! idListFileName.isEmpty());
QList<QString> idList;
if (isList)
{
if (! QFile::exists(idListFileName))
{
errorStr = "The ID list file does not exist: " + idListFileName;
errorStr = "The ID list does not exist: " + idListFileName;
return false;
}

Expand All @@ -3670,7 +3677,7 @@ bool Crit3DMeteoGridDbHandler::exportDailyDataCsv(QString &errorStr, bool isTPre

if (idList.size() == 0)
{
errorStr = "The ID list file is empty: " + idListFileName;
errorStr = "The ID list is empty: " + idListFileName;
return false;
}
}
Expand Down Expand Up @@ -3708,72 +3715,41 @@ bool Crit3DMeteoGridDbHandler::exportDailyDataCsv(QString &errorStr, bool isTPre
continue;
}

// header
// write header
QTextStream out(&outputFile);
out << "Date, Tmin (C), Tmax (C), Tavg (C), Prec (mm)";
if (isTPrec)
out << "\n";
else
out << "RHmin (%), RHmax (%), RHavg (%), Windspeed (m/s), Rad (MJ)\n";
out << "Date";
for (int i = 0; i < variableList.size(); i++)
{
if (variableList[i] != noMeteoVar)
{
std::string varName = getMeteoVarName(variableList[i]);
std::string unit = getUnitFromVariable(variableList[i]);
QString VarString = QString::fromStdString(varName + " (" + unit + ")");
out << "," + VarString;
}
}
out << "\n";

// save data
// write data
QDate currentDate = firstDate;
while (currentDate <= lastDate)
{
Crit3DDate myDate = getCrit3DDate(currentDate);
out << currentDate.toString("yyyy-MM-dd");

float tmin = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyAirTemperatureMin);
QString tminStr = "";
if (tmin != NODATA)
tminStr = QString::number(tmin);

float tmax = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyAirTemperatureMax);
QString tmaxStr = "";
if (tmax != NODATA)
tmaxStr = QString::number(tmax);

float tavg = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyAirTemperatureAvg);
QString tavgStr = "";
if (tavg != NODATA)
tavgStr = QString::number(tavg);

float prec = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyPrecipitation);
QString precStr = "";
if (prec != NODATA)
precStr = QString::number(prec);

out << currentDate.toString("yyyy-MM-dd") << "," << tminStr << "," << tmaxStr << "," << tavgStr << "," << precStr;
if (isTPrec)
out << "\n";
else
for (int i = 0; i < variableList.size(); i++)
{
float rhmin = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyAirRelHumidityMin);
QString rhminStr = "";
if (rhmin != NODATA)
rhminStr = QString::number(rhmin);

float rhmax = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyAirRelHumidityMax);
QString rhmaxStr = "";
if (rhmax != NODATA)
rhmaxStr = QString::number(rhmax);

float rhavg = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyAirRelHumidityAvg);
QString rhavgStr = "";
if (rhavg != NODATA)
rhavgStr = QString::number(rhavg);

float wspeed = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyWindScalarIntensityAvg);
QString wspeedStr = "";
if (wspeed != NODATA)
wspeedStr = QString::number(wspeed);

float rad = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, dailyGlobalRadiation);
QString radStr = "";
if (rad != NODATA)
radStr = QString::number(rad);

out << "," << rhminStr << "," << rhmaxStr << "," << rhavgStr << "," << wspeedStr << "," << radStr << "\n";
if (variableList[i] != noMeteoVar)
{
float value = _meteoGrid->meteoPointPointer(row,col)->getMeteoPointValueD(myDate, variableList[i]);
QString valueString = "";
if (value != NODATA)
valueString = QString::number(value);

out << "," << valueString;
}
}
out << "\n";

currentDate = currentDate.addDays(1);
}
Expand Down
Loading

0 comments on commit cd01594

Please sign in to comment.