/*
* This file is part of JGrasstools (http://www.jgrasstools.org)
* (C) HydroloGIS - www.hydrologis.com
*
* JGrasstools is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
package org.jgrasstools.hortonmachine.modules.networktools.trento_p.net;
import static java.lang.Math.abs;
import static java.lang.Math.pow;
import static java.lang.Math.sin;
import static java.lang.Math.sqrt;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.*;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.DEFAULT_CELERITY_FACTOR;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.DEFAULT_TMAX;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.EIGHTOVERTHREE;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.METER2CM;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.MINUTE2SEC;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.TWOOVERTHREE;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.TWO_THIRTEENOVERTHREE;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Utility.angleToFillDegree;
import java.text.DecimalFormat;
import java.text.NumberFormat;
import java.util.HashMap;
import java.util.LinkedHashMap;
import java.util.Map.Entry;
import java.util.Set;
import org.jgrasstools.gears.libs.modules.ModelsEngine;
import org.jgrasstools.gears.libs.monitor.IJGTProgressMonitor;
import org.jgrasstools.gears.utils.math.NumericsUtilities;
import org.jgrasstools.gears.utils.sorting.QuickSortAlgorithm;
import org.jgrasstools.hortonmachine.i18n.HortonMessageHandler;
import org.jgrasstools.hortonmachine.modules.networktools.trento_p.OmsTrentoP;
import org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Utility;
import org.joda.time.DateTime;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.OUT_ID_PIPE;
import static org.jgrasstools.hortonmachine.modules.networktools.trento_p.utils.Constants.OUT_INDEX_PIPE;
/**
*
* Class which contains methods to verify a sewer network.
*
* <p>
*
* It contains a geosewere method which verify the net..
*
* </p>
*
* @author Daniele Andreis, Riccardo Rigon, David Tamanini
*
*/
public class NetworkCalibration implements Network {
/**
* The initial time, it's the initial time whre to start to search the tp max.
* <p>
* It's fixed to 15 minutes.
*
* </p>
*
*/
public static final Integer INITIAL_TIME = 15;
DateTime first = null;
/*
* Time step.
*/
private Integer dt;
/*
* Monitor.
*/
private final IJGTProgressMonitor pm;
/*
*
*/
private final HortonMessageHandler msg = HortonMessageHandler.getInstance();
/*
* Rain data.
*/
private double[][] rainData;
/*
* The input rain.
*/
private final HashMap<DateTime, double[]> inRain;
/*
* Flag to exit if a pipe is full.
*/
boolean isFill = false;
/*
* Dati relativi alla rete
*/
private final Pipe[] networkPipes;
/*
* celerita-.
*/
private final double celerityfactor1;
/*
* The time of fill degree..
*/
private double[][] lastTimeDischarge;
/*
* The time of fill degree..
*/
private double[][] lastTimeFillDegree;
/*
* The end time of the simulation.
*/
private final int tMax;
/*
* A string where to put the warning messages.
*/
private final StringBuilder strBuilder;
/*
* The fill degree for each pipe and for each time.
*/
private HashMap<DateTime, HashMap<Integer, double[]>> fillDegree;
/*
* The discharg for each pipe and for each time.
*/
private HashMap<DateTime, HashMap<Integer, double[]>> discharge;
/*
* A flag that indicate to the program if search or not the maximum rain time.
*/
private final boolean foundMaxrainTime;
/*
* The upper limit time range. It is used to search the tpMax. The range to search the tpMax is [IITIAL_TIME, tpMaxCalibration].
*/
private final Integer tpMaxCalibration;
/*
* The rain time that give the maximum discharge.
*/
private int tpMax;
/*
* Is the number of data rain to use.
*/
private int nTime;
/*
* True if there is an infinite loop.
*/
private boolean infiniteLoop = false;
/*
* Max number of iteration to search the solution.
*/
private final static int MAX_NUMBER_ITERATION = 1000;
/**
* Builder for the Calibration class.
*/
public static class Builder {
// Mandatory parameters.
/*
* Progress monitor.
*/
private final IJGTProgressMonitor pm;
/*
* time step.
*/
private final Integer dt;
/*
* Builder that contains the warning messages.
*/
private final StringBuilder strBuilder;
/*
* The fillDegree for each pipe and for each time.
*/
private final HashMap<DateTime, HashMap<Integer, double[]>> fillDegree;
/*
* The discharg for each pipe and for each time.
*/
private final HashMap<DateTime, HashMap<Integer, double[]>> discharge;
/*
* Dati di pioggia.
*/
private final HashMap<DateTime, double[]> inRain;
/*
* Dati relativi alla rete
*/
private final Pipe[] networkPipe;
/*
* Flag to indicate if search or not the tpMax.
*/
private final boolean foundMaxrainTime;
/*
* Maximum time of calibration.
*/
private final Integer tpMaxCalibration;
// optional parameter.
/*
* Celerity.
*/
private double celerityfactor1 = DEFAULT_CELERITY_FACTOR;
/*
* max number of time step.
*/
private int tMax = (int) DEFAULT_TMAX;
/**
* Initialize the object with the needed parameters
*
*
* @param pm
* @param msg
* @param networkPipe the array of pipes.
* @param dt time step.
* @param inRain the rain data.
* @param outDischarge the output, discharge.
* @param outFillDegreethe output, fill degree.
* @param strBuilder a string used to store the warnings.
*/
public Builder( IJGTProgressMonitor pm, Pipe[] networkPipe, Integer dt, HashMap<DateTime, double[]> inRain,
HashMap<DateTime, HashMap<Integer, double[]>> outDischarge,
HashMap<DateTime, HashMap<Integer, double[]>> outFillDegree, StringBuilder strBuilder, Integer tpMaxCalibration,
boolean foundTpMax ) {
this.pm = pm;
this.networkPipe = networkPipe;
this.inRain = inRain;
this.dt = dt;
this.discharge = outDischarge;
this.fillDegree = outFillDegree;
this.strBuilder = strBuilder;
this.foundMaxrainTime = foundTpMax;
if (tpMaxCalibration != null) {
this.tpMaxCalibration = tpMaxCalibration;
} else {
this.tpMaxCalibration = tMax;
}
}
/**
* Set the max celerity factor.
*
* @param max
* celerityFactor.
*/
public Builder celerityFactor( double celerityFactor ) {
this.celerityfactor1 = celerityFactor;
return this;
}
/**
* Set the max number of time step.
*
* @param tMax
* .
*/
public Builder tMax( int tMax ) {
this.tMax = tMax;
return this;
}
/**
* Create a NetworkCalibration Object.
*
* @return
*/
public NetworkCalibration build() {
return new NetworkCalibration(this);
}
}
/**
* Set the parameter throughout the Builder,
*
*/
private NetworkCalibration( Builder builder ) {
this.dt = builder.dt;
this.pm = builder.pm;
this.celerityfactor1 = builder.celerityfactor1;
this.discharge = builder.discharge;
this.fillDegree = builder.fillDegree;
this.tMax = builder.tMax;
this.strBuilder = builder.strBuilder;
this.tpMaxCalibration = builder.tpMaxCalibration;
this.foundMaxrainTime = builder.foundMaxrainTime;
if (builder.networkPipe != null) {
this.networkPipes = builder.networkPipe;
} else {
pm.errorMessage(msg.message("trentoP.error.network"));
throw new IllegalArgumentException("trentoP.error.network");
}
if (builder.inRain != null) {
this.inRain = builder.inRain;
// create the rains array from the input.
Set<Entry<DateTime, double[]>> rainSet = inRain.entrySet();
DateTime second = null;
int l = rainSet.size();
rainData = new double[l][2];
int index = 0;
for( Entry<DateTime, double[]> rainRecord : rainSet ) {
DateTime dateTime = rainRecord.getKey();
double[] values = rainRecord.getValue();
if (first == null) {
first = dateTime;
} else if (second == null && dt == null) {
second = dateTime;
}
rainData[index][0] = index + 1;
rainData[index][1] = values[0];
index++;
}
// Evaluate the time step, if it isn't passed as a parameter,as a difference between two
// time.
if (dt == null) {
dt = abs(second.getMinuteOfDay() - first.getMinuteOfDay());
}
// if the input has the date in a wrong order.
if (dt <= 0) {
pm.errorMessage(msg.message("trentoP.error.t"));
throw new IllegalArgumentException(msg.message("trentoP.error.t"));
}
double tMaxApproximate = ModelsEngine.approximate2Multiple(tMax, dt);
// initialize the output.
nTime = (int) (tMaxApproximate / dt);
lastTimeDischarge = createMatrix();
lastTimeFillDegree = createMatrix();
} else {
pm.errorMessage(msg.message("trentoP.error.rainData "));
throw new IllegalArgumentException(msg.message("trentoP.error.rainData"));
}
}
/**
* Return the rain time.
*
* @return the rain time which give the maximum discharge.
*/
public int getTpMax() {
return tpMax;
}
/*
* Create the matrix which contains the result.
*/
private double[][] createMatrix() {
double time = 0;
double tmin = 0;
tmin = rainData[0][0];
time = tmin;
double[][] matrix = new double[nTime][networkPipes.length + 1];
for( int i = 0; i < nTime; ++i ) {
matrix[i][0] = time;
time += dt;
}
return matrix;
}
/**
* verify of the no-head pipes.
*
* <p>
* It evaluate the discharge.
* </p>
*
* @param k ID of the pipe where evaluate the discharge.
* @param cDelays delay matrix (for the evalutation of the flow wave).
* @param net matrix that contains value of the network.
* @return
*/
private double internalPipeVerify( int k, double[] cDelays, double[][] net, double[][] timeDischarge,
double[][] timeFillDegree, int tp ) {
int num;
double localdelay, olddelay, qMax, B, known, theta, u;
double[][] qPartial;
qPartial = new double[timeDischarge.length][timeDischarge[0].length];
calculateDelays(k, cDelays, net);
// First attempt local delay [min]
localdelay = 1;
double accuracy = networkPipes[0].getAccuracy();
int jMax = networkPipes[0].getjMax();
double minG = networkPipes[0].getMinG();
double maxtheta = networkPipes[0].getMaxTheta();
double tolerance = networkPipes[0].getTolerance();
int count = 0;
do {
olddelay = localdelay;
qMax = 0;
// Updates delays
for( int i = 0; i < net.length; i++ ) {
net[i][2] += localdelay;
};
for( int j = 0; j < net.length; ++j ) {
num = (int) net[j][0];
getHydrograph(num, qPartial, olddelay, net[j][2], tp);
}
getHydrograph(k, qPartial, olddelay, 0, tp);
qMax = ModelsEngine.sumDoublematrixColumns(k, qPartial, timeDischarge, 0, qPartial[0].length - 1, pm);
if (qMax <= 1)
qMax = 1;
// Resets delays
for( int i = 0; i < net.length; i++ ) {
net[i][2] -= localdelay;
}
calculateFillDegree(k, timeDischarge, timeFillDegree);
B = qMax / (CUBICMETER2LITER * networkPipes[k].getKs() * sqrt(networkPipes[k].verifyPipeSlope / METER2CM));
known = (B * TWO_THIRTEENOVERTHREE) / pow(networkPipes[k].diameterToVerify / METER2CM, EIGHTOVERTHREE);
theta = Utility.thisBisection(maxtheta, known, TWOOVERTHREE, minG, accuracy, jMax, pm, strBuilder);
// Average velocity in pipe [ m / s ]
u = qMax * 80 / (pow(networkPipes[k].diameterToVerify, 2) * (theta - sin(theta)));
localdelay = networkPipes[k].getLenght() / (celerityfactor1 * u * MINUTE2SEC);
count++;
// verify if it's an infiniteloop.
if (count > MAX_NUMBER_ITERATION) {
infiniteLoop = true;
throw new ArithmeticException();
}
} while( abs(localdelay - olddelay) / olddelay >= tolerance );
cDelays[k] = localdelay;
return qMax;
}
private void calculateFillDegree( int k, double[][] timeDischarge, double[][] timeFillDegree ) {
double accuracy = networkPipes[0].getAccuracy();
int jMax = networkPipes[0].getjMax();
double minG = networkPipes[0].getMinG();
double maxtheta = networkPipes[0].getMaxTheta();
double initialFillValue = angleToFillDegree(maxtheta) + 0.1;
for( int i = 0; i < timeDischarge.length; i++ ) {
// set it over the max value so if the bisection fails(because it's over the max) it's
// setted
timeFillDegree[i][k] = initialFillValue;
double q = timeDischarge[i][k];
if (q > NumericsUtilities.machineFEpsilon()) {
double B = q / (CUBICMETER2LITER * networkPipes[k].getKs() * sqrt(networkPipes[k].verifyPipeSlope / METER2CM));
double known = (B * TWO_THIRTEENOVERTHREE) / pow(networkPipes[k].diameterToVerify / METER2CM, EIGHTOVERTHREE);
double theta = Utility.thisBisection(maxtheta, known, TWOOVERTHREE, minG, accuracy, jMax, pm, strBuilder);
timeFillDegree[i][k] = angleToFillDegree(theta);
} else {
timeFillDegree[i][k] = 0.0;
}
}
}
/**
* Calcola il ritardo della tubazione k.
*
* @param k indice della tubazione.
* @param cDelays matrice dei ritardi.
* @param net matrice che contiene la sottorete.
*/
private void calculateDelays( int k, double[] cDelays, double[][] net )
{
double t;
int ind, r = 1;
for( int j = 0; j < net.length; ++j ) {
t = 0;
r = 1;
ind = (int) net[j][0];
/*
* Area k is not included in delays
*/
while( networkPipes[ind].getIndexPipeWhereDrain() != k ) {
ind = networkPipes[ind].getIndexPipeWhereDrain();
t += cDelays[ind];
r++;
}
if (r > networkPipes.length) {
pm.errorMessage(msg.message("trentoP.error.incorrectmatrix"));
throw new ArithmeticException(msg.message("trentoP.error.incorrectmatrix"));
}
net[j][2] = t;
}
}
/**
* Restituisce l'idrogramma.
*
* @param k tratto di tubazione.
* @param Qpartial matrice delle portate temporanee necessarie per costriure 'idrogramma.
* @param localdelay ritardo della tubazione k.
* @param delay ritardo temporale.
*/
private double getHydrograph( int k, double[][] Qpartial, double localdelay, double delay, int tp )
{
double Qmax = 0;
double tmin = rainData[0][0]; /* [min] */
int j = 0;
double t = tmin;
double Q;
double rain;
int maxRain = 0;
if (tMax == tpMaxCalibration) {
maxRain = rainData.length;
} else {
maxRain = tp / dt;
}
double tMaxApproximate = ModelsEngine.approximate2Multiple(tMax, dt);
for( t = tmin, j = 0; t <= tMaxApproximate; t += dt, ++j ) {
Q = 0;
for( int i = 0; i <= maxRain - 1; ++i ) {
// [ l / s ]
rain = rainData[i][1] * networkPipes[k].getDrainArea() * networkPipes[k].getRunoffCoefficient()
* HAOVERH_TO_METEROVERS;
if (t <= i * dt) {
Q += 0;
} else if (t <= (i + 1) * dt) {
Q += rain * pFunction(k, t - i * dt, localdelay, delay);
} else {
Q += rain * (pFunction(k, t - i * dt, localdelay, delay) - pFunction(k, t - (i + 1) * dt, localdelay, delay));
}
}
Qpartial[j][k] = Q;
if (Q >= Qmax) {
Qmax = Q;
}
}
return Qmax;
}
/**
*
* Verify of the head pipes.
*
* <p>
* It evaluate the discharge.
* </p>
*
* @param k
* ID of the pipe where evaluate the discharge.
* @param cDelays
* delay matrix (for the evalutation of the flow wave).
*/
private double headPipeVerify( int k, double[] cDelays, double[][] timeDischarge, double[][] timeFillDegree, int tp ) {
double olddelay = 0;
double qMax = 0;
double B = 0;
double known = 0;
double theta = 0;
double u = 0;
/* First attempt local delay [min] */
double localdelay = 1;
double accuracy = networkPipes[0].getAccuracy();
int jMax = networkPipes[0].getjMax();
double minG = networkPipes[0].getMinG();
double maxtheta = networkPipes[0].getMaxTheta();
double tolerance = networkPipes[0].getTolerance();
int count = 0;
do {
olddelay = localdelay;
qMax = getHydrograph(k, timeDischarge, olddelay, 0, tp);
if (qMax <= 1) {
qMax = 1;
}
calculateFillDegree(k, timeDischarge, timeFillDegree);
B = qMax / (CUBICMETER2LITER * networkPipes[k].getKs() * Math.sqrt(networkPipes[k].verifyPipeSlope / METER2CM));
known = (B * TWO_THIRTEENOVERTHREE) / Math.pow(networkPipes[k].diameterToVerify / METER2CM, EIGHTOVERTHREE);
theta = Utility.thisBisection(maxtheta, known, TWOOVERTHREE, minG, accuracy, jMax, pm, strBuilder);
double tmp1 = 0;
double tmp2 = 0;
if (k - 1 >= 0) {
tmp1 = networkPipes[k].diameterToVerify;
tmp2 = networkPipes[k].getLenght();
}
// Average velocity in pipe [ m / s ]
u = qMax * 80 / (Math.pow(tmp1, 2) * (theta - Math.sin(theta)));
localdelay = tmp2 / (celerityfactor1 * u * MINUTE2SEC);
count++;
if (count > MAX_NUMBER_ITERATION) {
infiniteLoop = true;
throw new ArithmeticException();
}
} while( Math.abs(localdelay - olddelay) / olddelay >= tolerance );
cDelays[k] = localdelay;
return qMax;
}
/**
*
* restituisce la funzione p.
*
* @param k tratto di tubazione esaminata.
* @param ttempo esaminato.
* @param localdelay ritardo nella tubazione k.
* @param delay ritardo totale.
* @return il valore della funxione p
*/
private double pFunction( int k, double t, double localdelay, double delay ) {
double P = 0;
if (t < 0) {
pm.errorMessage(msg.message("trentoP.error.negativeP"));
throw new ArithmeticException(msg.message("trentoP.error.negativeP"));
}
if (t < delay) {
P = 0;
} else if (t <= (delay + localdelay)) {
P = (t - delay) / localdelay + networkPipes[k].k / localdelay * (Math.exp(-(t - delay) / networkPipes[k].k) - 1);
} else {
P = 1 + networkPipes[k].k / localdelay * Math.exp(-(t - delay) / networkPipes[k].k)
* (1 - Math.exp(localdelay / networkPipes[k].k));
}
return P;
}
/**
* Estimate the discharge for each time and for each pipes.
*
*<p>
*It can work with a single rain time step (if there is an actual rain) or search the maximum rain time, if the rain is unknown and
*the rain data are builtthroghout the rain possibility curve.
*</p>
* <p>
* It work throgout 2 loop:
* <ol>
* *
* <li>The first on the head pipes.
* <li>The second on the pipes which are internal.
*
* </ol>
* </p>
*/
@Override
public void geoSewer() throws Exception {
if (!foundMaxrainTime) {
evaluateDischarge(lastTimeDischarge, lastTimeFillDegree, tpMaxCalibration);
} else {
/*
* start to evaluate the discharge from 15 minutes,evaluate the nearsted value to 15 minutes.
*/
int minTime = (int) ModelsEngine.approximate2Multiple(INITIAL_TIME, dt);
double qMax = 0;
for( int i = minTime; i < tpMaxCalibration; i = i + dt ) {
tpMax = i;
double[][] timeDischarge = createMatrix();
double[][] timeFillDegree = createMatrix();
double q = evaluateDischarge(timeDischarge, timeFillDegree, i);
if (q > qMax) {
qMax = q;
lastTimeDischarge = timeDischarge;
lastTimeFillDegree = timeFillDegree;
} else if (q < qMax) {
break;
}
if (isFill) {
break;
}
}
}
getNetData();
}
/*
* Fill the two output HashMap.
*/
private void getNetData() {
int nTime = lastTimeDischarge.length;
int length = lastTimeDischarge[0].length;
HashMap<Integer, double[]> tmpHMDis = new LinkedHashMap<Integer, double[]>();
HashMap<Integer, double[]> tmpHMFill = new LinkedHashMap<Integer, double[]>();
// order the outpt.
int netLength = networkPipes.length;
double[] one = new double[netLength];
double[] two = new double[netLength];
for( int i = 0; i < netLength; i++ ) {
one[i] = i;
two[i] = networkPipes[i].getId();
}
QuickSortAlgorithm sort = new QuickSortAlgorithm(pm);
sort.sort(two, one);
for( int i = 0; i < length - 1; i++ ) {
int index = (int) one[i];
tmpHMDis.put(networkPipes[index].getId(), new double[]{lastTimeDischarge[0][index]});
tmpHMFill.put(networkPipes[index].getId(), new double[]{lastTimeFillDegree[0][index]});
}
discharge.put(first, tmpHMDis);
fillDegree.put(first, tmpHMFill);
DateTime tmp = first;
for( int i = 1; i < nTime; ++i ) {
tmp = tmp.plusMinutes(dt);
tmpHMDis = new LinkedHashMap<Integer, double[]>();
tmpHMFill = new LinkedHashMap<Integer, double[]>();
for( int j = 0; j < length - 1; j++ ) {
int index = (int) one[j];
tmpHMDis.put(networkPipes[index].getId(), new double[]{lastTimeDischarge[i][index]});
tmpHMFill.put(networkPipes[index].getId(), new double[]{lastTimeFillDegree[i][index]});
}
discharge.put(tmp, tmpHMDis);
fillDegree.put(tmp, tmpHMFill);
}
}
/**
* Compila la mantrice net con tutte i dati del sottobacino con chiusura nel
* tratto che si sta analizzando, e restituisce la sua superfice
*
* @param k tratto analizzato del sottobacino chiuso in l.
* @param l chiusura del bacino.
* @param one indice dei versanti.
* @param net sottobacino che si chiude in l.
*/
private double scanNetwork( int k, int l, double[] one, double[][] net ) {
int ind;
/*
* t Ritardo accumulato dall'onda prima di raggiungere il tratto si sta
* dimensionando.
*/
double t;
/*
* Distanza percorsa dall'acqua dall'area dove e' caduta per raggiungere
* l'ingresso del tratto che si sta dimensionando.
*/
double length;
/*
* Superfice del sottobacino con chiusura nel tratto che si sta
* analizzando.
*/
double totalarea = 0;
int r = 0;
int i = 0;
/*
* In one gli stati sono in ordine di magmitude crescente. Per ogni
* stato di magnitude inferiore a quella del tratto l che si sta
* progettando.
*/
for( int j = 0; j < k; j++ ) {
/* La portata e valutata all'uscita di ciascun tratto */
t = 0;
/*
* ID dello lo stato di magnitude inferiore a quello del tratto che
* si sta progettando.
*/
i = (int) one[j];
ind = i;
// la lunghezza del tubo precedentemente progettato
length = networkPipes[ind].getLenght();
// seguo il percorso dell'acqua finch� non si incontra l'uscita.
while( networkPipes[ind].getIdPipeWhereDrain() != OUT_ID_PIPE ) {
// lo stato dove drena a sua volta.
ind = networkPipes[ind].getIndexPipeWhereDrain();
/*
* se lo stato drena direttamente in quello che si sta
* progettando
*/
if (ind == l) {
/*
* ID dello stato che drena in l, piu o meno direttamente.
*/
net[r][0] = i;
/*
* lunghezza del percorsa dall'acqua prima di raggiungere lo
* stato l che si sta progettando
*/
net[r][1] = length + networkPipes[l].getLenght();
/*
* Ritardo accumulato dall'onda di piena formatasi in uno
* degli stati a monte, prima di raggiungere il tratto l che
* si sta progettando
*/
net[r][2] = t;
/*
* area di tutti gli stati a monte che direttamente o
* indirettamente drenano in l
*/
totalarea += networkPipes[i].getDrainArea();
r++;
break;
}
}
/*
* viene incrementato solo se l'area drena in l quindi non puo'
* superare net->nrh
*/
if (r > net.length)
break;
}
// area degli stati a monte che drenano in l, l compreso
totalarea += networkPipes[l].getDrainArea();
return totalarea;
}
private double evaluateDischarge( double[][] timeDischarge, double[][] timeFillDegree, int tp ) {
/* l Tratto che si sta progettando. */
int l;
/*
* matrice che per ciascun area non di testa, contiene i dati geometrici
* degli stati a monte, che direttamente o indirettamente, drenano in
* esso
*/
double[][] net;
/*
* contiene la magnitude dei vari stati.
*/
double[] magnitude = new double[networkPipes.length];
/*
* vettore che contiene l'indice dei versanti.
*/
double[] one = new double[networkPipes.length];
/*
* vettore che contiene gli stati riceventi, compresa almeno un'uscita
*/
double[] two = new double[networkPipes.length];
/*
*Max discharge at this iteration;
*/
double qMax = 0;
// initialize the discharge array
for( int i = 0; i < networkPipes.length; i++ ) {
/* Indice degli stati */
one[i] = i;
/* Indice degli stati riceventi, compresa almeno un'uscita */
two[i] = networkPipes[i].getIndexPipeWhereDrain();
}
/* Calcola la magnitude di ciascun stato */
Utility.pipeMagnitude(magnitude, two, pm);/*
* Calcola la magnitude di
* ciascun stato
*/
/* al vettore two vengono assegnati gli elementi di magnitude */
for( int i = 0; i < two.length; i++ ) /*
* al vettore two vengono assegnati
* gli elementi di magnitude
*/
{
two[i] = magnitude[i];
}
/*
* Ordina gli elementi del vettore magnitude in ordine crescente, e
* posiziona nello stesso ordine gli elementi di one
*/
QuickSortAlgorithm t = new QuickSortAlgorithm(pm);
t.sort(magnitude, one);
int k = 0;
// tratto che si sta analizzando o progettando
l = (int) one[k];
pm.beginTask(msg.message("trentoP.begin"), networkPipes.length - 1);
isFill = false;
double[] cDelays = new double[networkPipes.length];
double maxFill = angleToFillDegree(networkPipes[k].getMaxTheta());
while( magnitude[k] == 1 ) {
try {
double q = headPipeVerify(l, cDelays, timeDischarge, timeFillDegree, tp);
if (q > qMax) {
qMax = q;
}
// Passo allo stato successivo
k++;
if (k < magnitude.length) {
/*
* Il prossimo tratto da progettare, ovviamente se avra
* magnitude=1
*/
l = (int) one[k];
} else {
break;
}
pm.worked(1);
} catch (ArithmeticException e) {
if (infiniteLoop) {
strBuilder.append(msg.message("trentoP.error.infiniteLoop"));
} else {
NumberFormat formatter = new DecimalFormat("#.###");
String limit = formatter.format(maxFill);
strBuilder.append(" ");
strBuilder.append(msg.message("trentoP.warning.emptydegree")); //$NON-NLS-2$
strBuilder.append(limit);
strBuilder.append(" ");
strBuilder.append(msg.message("trentoP.warning.emptydegree2"));
strBuilder.append(networkPipes[l - 1].getId());
strBuilder.append(" ");
strBuilder.append("tp " + tp);
strBuilder.append("\n");
isFill = true;
}
break;
}
}
/*
* ----- INIZIO CICLO WHILE PER LA PROGETTAZIONE DELLE AREE NON DI TESTA
* -----
*
* Magnitude > 1 AREE NON DI TESTA
*/
if (!isFill) {
while( k < magnitude.length ) {
try {
net = new double[(int) (magnitude[k] - 1)][9];
scanNetwork(k, l, one, net);
double q = internalPipeVerify(l, cDelays, net, timeDischarge, timeFillDegree, tp);
if (q > qMax) {
qMax = q;
}
/* Passo allo stato successivo */
k++;
/* se non sono arrivato alla fine */
if (k < magnitude.length) {
/* Prossimo stato da progettare */
l = (int) one[k];
} else {
break;
}
pm.worked(1);
} catch (ArithmeticException e) {
// if there is an infinite loop.
if (infiniteLoop) {
strBuilder.append(msg.message("trentoP.error.infiniteLoop"));
} else {
// if a pipe is fill.
strBuilder.append(msg.message("trentoP.warning.emptydegree")); //$NON-NLS-2$
strBuilder.append(maxFill);
strBuilder.append(" ");
strBuilder.append(msg.message("trentoP.warning.emptydegree2"));
strBuilder.append(networkPipes[l].getId());
strBuilder.append(" ");
strBuilder.append("tp " + tp);
strBuilder.append("\n");
isFill = true;
}
break;
}
}
}
return qMax;
}
}