Examples with EIAltimetry org.hortonmachine.gears.io.eicalculator.EIAltimetry used on opensource projects

Example 1 with EIAltimetry

use of org.hortonmachine.gears.io.eicalculator.EIAltimetry in project hortonmachine by TheHortonMachine.

the class OmsJami method process.

@Execute
public void process() throws Exception {
    pm.message("processing " + tCurrent + " " + pType);
    checkNull(inAltimetry, inAreas, inMeteo, inStations);
    // check fascie num
    int fascieNum = 0;
    for (EIAreas area : inAreas) {
        int idAltimetricBand = area.altimetricBandId;
        if (idAltimetricBand > fascieNum) {
            fascieNum = idAltimetricBand;
        }
    }
    // increment by one. Fascie start from 0, so 0-5 are 6 fascie
    fascieNum++;
    currentTimestamp = formatter.parseDateTime(tCurrent);
    outInterpolatedBand = new HashMap<Integer, double[]>();
    outInterpolated = new HashMap<Integer, double[]>();
    /*
         * get stations
         */
    pm.message("Read stations data.");
    stationCoordinates = new ArrayList<Coordinate>();
    stationFeatures = new ArrayList<SimpleFeature>();
    FeatureIterator<SimpleFeature> featureIterator = inStations.features();
    while (featureIterator.hasNext()) {
        SimpleFeature feature = featureIterator.next();
        Coordinate stationCoord = ((Geometry) feature.getDefaultGeometry()).getCoordinate();
        stationCoordinates.add(stationCoord);
        stationFeatures.add(feature);
    }
    featureIterator.close();
    pm.done();
    /*
         * get basins and create of every basin a buffered geometry at 10,
         * 20 and 50 km, which will be used to find the nearest stations
         * around.
         */
    basinBaricenterCoordinates = new ArrayList<Coordinate>();
    basinFeatures = new ArrayList<SimpleFeature>();
    featureIterator = inInterpolate.features();
    pm.beginTask("Read basins data.", inInterpolate.size());
    while (featureIterator.hasNext()) {
        pm.worked(1);
        SimpleFeature feature = featureIterator.next();
        Coordinate baricenterCoord = ((Geometry) feature.getDefaultGeometry()).getCentroid().getCoordinate();
        basinBaricenterCoordinates.add(baricenterCoord);
        basinFeatures.add(feature);
    }
    featureIterator.close();
    pm.done();
    statElev = new double[stationCoordinates.size()];
    statId = new double[stationCoordinates.size()];
    stationId2CoordinateMap = new HashMap<Integer, Coordinate>();
    extractFromStationFeatures();
    stationid2StationindexMap = new HashMap<Integer, Integer>();
    for (int i = 0; i < statId.length; i++) {
        stationid2StationindexMap.put((int) statId[i], i);
    }
    /*
         * get the basin's id attribute index
         */
    SimpleFeature tmpfeature = basinFeatures.get(0);
    SimpleFeatureType featureType = tmpfeature.getFeatureType();
    basinIdFieldIndex = featureType.indexOf(fBasinid);
    if (basinIdFieldIndex == -1) {
        throw new IllegalArgumentException("The field of the basin id couldn't be found in the supplied basin data.");
    }
    basinid2BasinindexMap = new HashMap<Integer, Integer>();
    basinindex2basinidMap = new HashMap<Integer, Integer>();
    for (int i = 0; i < basinBaricenterCoordinates.size(); i++) {
        int basinid = ((Number) basinFeatures.get(i).getAttribute(basinIdFieldIndex)).intValue();
        basinid2BasinindexMap.put(basinid, i);
        basinindex2basinidMap.put(i, basinid);
    }
    calculateAreas(fascieNum);
    pm.message("Creating the band's elevation for every basin matrix.");
    /*
         * create the altimetric bands matrix
         */
    int basinNum = basinBaricenterCoordinates.size();
    bandsBasins = new double[fascieNum][basinNum];
    for (int i = 0; i < inAltimetry.size(); i++) {
        EIAltimetry eiAltimetry = inAltimetry.get(i);
        int idbasin = eiAltimetry.basinId;
        int idfascia = eiAltimetry.altimetricBandId;
        double elevationInBandBaricenter = eiAltimetry.elevationValue;
        Integer index = basinid2BasinindexMap.get(idbasin);
        if (index != null)
            bandsBasins[idfascia][index] = elevationInBandBaricenter;
    // TODO make it range aware
    // double range = altimetryScalarSet.get(i + 3);
    // bandsBasins[idfascia +
    // 1][basinid2BasinindexMap.get(idbasin)] = baricenter
    // + range / 2.0;
    }
    double[] stationBinsArrays = new double[pBins + 1];
    double maxStatElev = statElev[statElev.length - 1];
    double minStatElev = statElev[0];
    double interval = (maxStatElev - minStatElev) / (double) pBins;
    for (int i = 0; i < stationBinsArrays.length; i++) {
        stationBinsArrays[i] = minStatElev + i * interval;
    }
    /*
         * find all stations inside a elevation band for every basin
         */
    pm.beginTask("Finding all stations inside a elevation band for every basin.", stationBinsArrays.length - 1);
    bin2StationsListMap = new HashMap<Integer, List<Integer>>();
    for (int i = 0; i < stationBinsArrays.length - 1; i++) {
        List<Integer> stationsIds = new ArrayList<Integer>();
        for (int j = 0; j < statId.length; j++) {
            double id = statId[j];
            double elev = statElev[j];
            if (elev >= stationBinsArrays[i] && elev < stationBinsArrays[i + 1]) {
                stationsIds.add((int) id);
            }
        }
        bin2StationsListMap.put(i, stationsIds);
        pm.worked(1);
    }
    pm.done();
    /*
         * get values for current timestep and order them with the stations ids
         */
    double[] statValues = new double[stationCoordinates.size()];
    for (int i = 0; i < statValues.length; i++) {
        statValues[i] = doubleNovalue;
    }
    Set<Integer> stationIdSet = inMeteo.keySet();
    for (Integer stationId : stationIdSet) {
        int id = stationId;
        double[] value = inMeteo.get(id);
        Integer index = stationid2StationindexMap.get((int) id);
        if (index == null)
            continue;
        statValues[index] = value[0];
    }
    // number of active stations for every basin
    int[] activeStationsPerBasin = new int[basinBaricenterCoordinates.size()];
    int[][] stazBacMatrix = createStationBasinsMatrix(statValues, activeStationsPerBasin);
    int[][] stations = new int[stazBacMatrix.length][stazBacMatrix[0].length];
    int contStations = 0;
    // riempimento vettori/matrici
    for (int i = 0; i < stazBacMatrix[0].length; i++) {
        // indice bacino
        contStations = 0;
        for (int j = 0; j < stazBacMatrix.length; j++) {
            // indice stazione
            if (stazBacMatrix[j][i] == 1) {
                stations[contStations][i] = j;
                contStations += 1;
            }
        }
    }
    int bandsNum = bandsBasins.length;
    if (pType == DTDAY || pType == DTMONTH) {
        /*
             * calculate the DT month and day for each station
             */
        // System.out.println("Calculating the dayly and monthly Dt for each station...");
        rangeT(statValues);
    }
    pm.beginTask("Interpolating over bands and basins...", basinBaricenterCoordinates.size());
    for (int i = 0; i < basinBaricenterCoordinates.size(); i++) {
        pm.worked(1);
        // interpolated value for every band
        double[] interpolatedMeteoForBand = new double[bandsNum];
        double interpolatedMeteoForBasin = 0;
        int cont = 0;
        double h;
        int[] jj_av;
        // trova le stazioni che forniscono dati
        // costruisco un nuovo
        jj_av = new int[activeStationsPerBasin[i]];
        // stazioni del bacino in studio
        for (int j = 0; j < activeStationsPerBasin[i]; j++) {
            if (pType != DTDAY && pType != DTMONTH) {
                if (!isNovalue(statValues[stations[j][i]])) {
                    // registro le stazioni
                    jj_av[cont] = stations[j][i];
                    // attive
                    cont += 1;
                }
            } else {
                // NODATA
                if (!isNovalue(minTempPerStation[stations[j][i]]) && isNovalue(maxTempPerStation[stations[j][i]])) {
                    // jj conterrà le stazioni che hanno dati di escursione
                    // termica
                    // giornaliera
                    // registro le stazioni
                    jj_av[cont] = stations[j][i];
                    // attive
                    cont += 1;
                }
            }
        }
        // sopravviva
        if (cont == 0) {
            if (pType == TEMPERATURE) {
                // caso dei dati di temperatura
                pm.errorMessage("ERRORE: PER IL BACINO " + i + " NON SONO DISPONIBILI DATI DI TEMPERATURA, PER QUESTO BACINO STAND-BY");
                for (int f = 0; f < bandsNum; f++) {
                    // per tutte le fasce
                    // altimetriche metto il
                    // dato a -100
                    interpolatedMeteoForBand[f] = doubleNovalue;
                }
                interpolatedMeteoForBasin = doubleNovalue;
            } else if (pType == PRESSURE) {
                // caso dei dati di pressione
                pm.message("  -> Per il bacino " + i + " non sono disponibili dati di pressione, uso valori di default");
                for (int f = 0; f < bandsNum; f++) {
                    // per tutte le fasce
                    // altimetriche considero
                    // un'adiabatica
                    interpolatedMeteoForBand[f] = 1013.25 * Math.exp(-(bandsBasins[f][i]) * 0.00013);
                    interpolatedMeteoForBasin = interpolatedMeteoForBasin + interpolatedMeteoForBand[f] * basinAreasPerFascias[i][f] / basinAreas[i];
                }
            } else if (pType == HUMIDITY) {
                // caso dei dati di umidità
                pm.message("  -> Per il bacino " + i + " non sono disponibili dati di umidita', uso valori di default");
                for (int f = 0; f < bandsNum; f++) {
                    // per tutte le fasce
                    // altimetriche metto NODATA
                    interpolatedMeteoForBand[f] = defaultRh;
                }
                interpolatedMeteoForBasin = defaultRh;
            } else if (pType == WIND) {
                // caso dei dati di velocità del vento
                pm.message("  -> Per il bacino " + i + " non sono disponibili dati di velocita' del vento, uso valori di default");
                for (int f = 0; f < bandsNum; f++) {
                    // per tutte le fasce
                    // altimetriche metto NODATA
                    interpolatedMeteoForBand[f] = defaultW;
                }
                interpolatedMeteoForBasin = defaultW;
            } else if (pType == DTDAY) {
                // caso dei dati di escursione termica
                // giornaliera
                pm.message("  -> Per il bacino " + i + " non sono disponibili dati di escursione termica giornaliera', uso valori di default");
                for (int f = 0; f < bandsNum; f++) {
                    // per tutte le fasce
                    // altimetriche del bacino
                    // assegno all'escursione termica giornaliera il dato
                    // DTd
                    // messo nel file dei parametri
                    interpolatedMeteoForBand[f] = defaultDtday;
                }
                interpolatedMeteoForBasin = defaultDtday;
            } else if (pType == DTMONTH) {
                // caso dei dati di escursione termica
                // mensile
                pm.message("  -> Per il bacino " + i + " non sono disponibili dati di escursione termica mensile', uso valori di default");
                for (int f = 0; f < bandsNum; f++) {
                    /*
                         *  per tutte le fasce
                         * altimetriche del bacino
                         */
                    // assegno all'escursione termica media mensile il
                    // datoDTm
                    // messo nel file dei parametri
                    interpolatedMeteoForBand[f] = defaultDtmonth;
                }
                interpolatedMeteoForBasin = defaultDtmonth;
            }
        } else if (cont == 1) {
            // standard per T e P, valori costanti per RH e V
            for (int f = 0; f < bandsNum; f++) {
                // altimetriche
                if (pType == TEMPERATURE) {
                    // trasformo la temp in K e calcolo T
                    // con
                    // l'adiabatica semplice
                    interpolatedMeteoForBand[f] = (statValues[jj_av[0]] + tk) * Math.exp(-(bandsBasins[f][i] - statElev[jj_av[0]]) * GAMMA / (statValues[jj_av[0]] + tk)) - tk;
                } else if (pType == PRESSURE) {
                    // calcolo P con il gradiente
                    // adiabatico
                    interpolatedMeteoForBand[f] = statValues[jj_av[0]] * Math.exp(-(bandsBasins[f][i] - statElev[jj_av[0]]) * 0.00013);
                } else if (pType == DTDAY) {
                    // se ho una sola stazione assegno il valore della
                    // stazione a tutto il
                    // bacino
                    // altimetriche del bacino assegno il valore di
                    // escursione massima
                    // giornaliera
                    interpolatedMeteoForBand[f] = maxTempPerStation[jj_av[0]] - minTempPerStation[jj_av[0]];
                    if ((maxTempPerStation[jj_av[0]] - minTempPerStation[jj_av[0]]) <= 0) {
                        interpolatedMeteoForBand[f] = defaultDtday;
                    }
                } else if (pType == DTMONTH) {
                    // se ho una sola stazione assegno il valore della
                    // stazione a tutto il
                    // bacino
                    // altimetriche del bacino assegno il valore di
                    // escursione massima mensile
                    interpolatedMeteoForBand[f] = DTmonth[jj_av[0]];
                } else {
                    // RH e V sono costanti al variare delle fasce
                    // altimetriche
                    interpolatedMeteoForBand[f] = statValues[jj_av[0]];
                }
                interpolatedMeteoForBasin = interpolatedMeteoForBasin + interpolatedMeteoForBand[f] * basinAreasPerFascias[i][f] / basinAreas[i];
            }
        } else {
            // caso 2. ci sono almeno 2 stazioni (a quote inferiori alla
            // stazioni piu' bassa considero atmosfera standard come a quote
            // superiori alla staz. piu' alta, in mezzo calcolo LAPSE RATE)
            // alloca L (vettore di dimensioni numero di stazioni attive-1)
            double[] lapseRate = new double[cont - 1];
            for (int j = 0; j < cont - 1; j++) {
                // le stazioni sono in
                // ordine di
                // quota
                // L[j] e' il lapse rate tra la stazione j e j+1, puo'
                // essere
                // calcolato dai dati per j che va da 1 a n-1, dove n e' il
                // numero di stazioni (cont)
                lapseRate[j] = (statValues[jj_av[j]] - statValues[jj_av[j + 1]]) / (statElev[jj_av[j + 1]] - statElev[jj_av[j]]);
            }
            for (int f = 0; f < bandsNum; f++) {
                // più bassa
                if (bandsBasins[f][i] <= statElev[jj_av[0]]) {
                    if (pType == TEMPERATURE) {
                        // T
                        interpolatedMeteoForBand[f] = statValues[jj_av[0]] - GAMMA * (bandsBasins[f][i] - statElev[jj_av[0]]);
                    } else if (pType == PRESSURE) {
                        // P
                        interpolatedMeteoForBand[f] = statValues[jj_av[0]] - (statValues[jj_av[0]] * 0.00013) * (bandsBasins[f][i] - statElev[jj_av[0]]);
                    } else if (pType == DTDAY) {
                        interpolatedMeteoForBand[f] = maxTempPerStation[jj_av[0]] - minTempPerStation[jj_av[0]];
                        if ((maxTempPerStation[jj_av[0]] - minTempPerStation[jj_av[0]]) <= 0) {
                            interpolatedMeteoForBand[f] = defaultDtday;
                        }
                    } else if (pType == DTMONTH) {
                        interpolatedMeteoForBand[f] = DTmonth[jj_av[0]];
                    } else {
                        // RH e V
                        interpolatedMeteoForBand[f] = statValues[jj_av[0]];
                    }
                // per le fasce altimetriche con quote piu' alte della
                // quota
                // della stazione piu' alta prendo i dati della stazione
                // più alta
                } else if (bandsBasins[f][i] >= statElev[jj_av[cont - 1]]) {
                    if (pType == TEMPERATURE) {
                        // T
                        interpolatedMeteoForBand[f] = statValues[jj_av[cont - 1]] - GAMMA * (bandsBasins[f][i] - statElev[jj_av[cont - 1]]);
                    } else if (pType == PRESSURE) {
                        // P
                        interpolatedMeteoForBand[f] = statValues[jj_av[cont - 1]] - (statValues[jj_av[cont - 1]] * 0.00013) * (bandsBasins[f][i] - statElev[jj_av[cont - 1]]);
                    } else if (pType == DTDAY) {
                        interpolatedMeteoForBand[f] = maxTempPerStation[jj_av[cont - 1]] - minTempPerStation[jj_av[cont - 1]];
                        if ((maxTempPerStation[jj_av[0]] - minTempPerStation[jj_av[0]]) <= 0) {
                            interpolatedMeteoForBand[f] = defaultDtday;
                        }
                    } else if (pType == DTMONTH) {
                        interpolatedMeteoForBand[f] = DTmonth[jj_av[cont - 1]];
                    } else {
                        // RH e V
                        interpolatedMeteoForBand[f] = statValues[jj_av[cont - 1]];
                    }
                } else {
                    int k = cont - 1;
                    if (pType == DTDAY) {
                        // la max delle stazioni
                        do {
                            k -= 1;
                            h = statElev[jj_av[k]];
                        } while (bandsBasins[f][i] <= h);
                        // for (int j = 0; j < cont; j++) {
                        // if (f ==0 && i == 100) {
                        // System.out.println(j + " "+ statElev[jj_av[j]]);
                        // }
                        // }
                        // interpolatedMeteoForBand[f] =
                        // ((maxTempPerStation[jj_av[k]] -
                        // minTempPerStation[jj_av[k]])
                        // * (statElev[jj_av[k + 1]] - bandsBasins[f][i]) +
                        // (maxTempPerStation[jj_av[k + 1]] -
                        // minTempPerStation[jj_av[k + 1]])
                        // * (bandsBasins[f][i] - statElev[jj_av[k]]))
                        // / (statElev[jj_av[k + 1]] - statElev[jj_av[k]]);
                        interpolatedMeteoForBand[f] = ((maxTempPerStation[jj_av[k + 1]] - minTempPerStation[jj_av[k + 1]]) - (maxTempPerStation[jj_av[k]] - minTempPerStation[jj_av[k]])) * (bandsBasins[f][i] - statElev[jj_av[k]]) / (statElev[jj_av[k + 1]] - statElev[jj_av[k]]) + (maxTempPerStation[jj_av[k]] - minTempPerStation[jj_av[k]]);
                        if (interpolatedMeteoForBand[f] <= 0) {
                            interpolatedMeteoForBand[f] = defaultDtday;
                        }
                    } else if (pType == DTMONTH) {
                        // la max delle stazioni
                        do {
                            k -= 1;
                            h = statElev[jj_av[k]];
                        } while (bandsBasins[f][i] <= h);
                        interpolatedMeteoForBand[f] = (DTmonth[jj_av[k]] * (statElev[jj_av[k + 1]] - bandsBasins[f][i]) + DTmonth[jj_av[k + 1]] * (bandsBasins[f][i] - statElev[jj_av[k]])) / (statElev[jj_av[k + 1]] - statElev[jj_av[k]]);
                    } else {
                        do {
                            k -= 1;
                            h = statElev[jj_av[k]];
                        } while (bandsBasins[f][i] <= h);
                        interpolatedMeteoForBand[f] = statValues[jj_av[k]] - lapseRate[k] * (bandsBasins[f][i] - statElev[jj_av[k]]);
                    }
                }
                interpolatedMeteoForBasin = interpolatedMeteoForBasin + interpolatedMeteoForBand[f] * basinAreasPerFascias[i][f] / basinAreas[i];
            }
            // controllo su RH>100 e v=0
            if (pType == HUMIDITY) {
                // RH
                double MAX_HUMIDITY = 100;
                double MIN_HUMIDITY = 5;
                for (int f = 0; f < bandsNum; f++) {
                    if (interpolatedMeteoForBand[f] > MAX_HUMIDITY)
                        interpolatedMeteoForBand[f] = MAX_HUMIDITY;
                    if (interpolatedMeteoForBand[f] < MIN_HUMIDITY)
                        interpolatedMeteoForBand[f] = MIN_HUMIDITY;
                }
                if (interpolatedMeteoForBasin > MAX_HUMIDITY)
                    interpolatedMeteoForBasin = MAX_HUMIDITY;
                if (interpolatedMeteoForBasin < MIN_HUMIDITY)
                    interpolatedMeteoForBasin = MIN_HUMIDITY;
            } else if (pType == WIND) {
                // V
                double MIN_WIND = 0.01;
                for (int f = 0; f < bandsNum; f++) {
                    if (interpolatedMeteoForBand[f] < MIN_WIND)
                        interpolatedMeteoForBand[f] = MIN_WIND;
                }
                if (interpolatedMeteoForBasin < MIN_WIND)
                    interpolatedMeteoForBasin = MIN_WIND;
            }
        }
        int basinid = ((Number) basinFeatures.get(i).getAttribute(basinIdFieldIndex)).intValue();
        outInterpolatedBand.put(basinid, interpolatedMeteoForBand);
        outInterpolated.put(basinid, new double[] { interpolatedMeteoForBasin });
    }
    pm.done();
}

Example 2 with EIAltimetry

use of org.hortonmachine.gears.io.eicalculator.EIAltimetry in project hortonmachine by TheHortonMachine.

the class OmsEnergyIndexCalculator method process.

@Execute
public void process() throws Exception {
    HashMap<String, Double> regionMap = CoverageUtilities.getRegionParamsFromGridCoverage(inBasins);
    cols = regionMap.get(CoverageUtilities.COLS).intValue();
    rows = regionMap.get(CoverageUtilities.ROWS).intValue();
    dx = regionMap.get(CoverageUtilities.XRES);
    dy = regionMap.get(CoverageUtilities.YRES);
    double n = regionMap.get(CoverageUtilities.NORTH);
    double s = regionMap.get(CoverageUtilities.SOUTH);
    double w = regionMap.get(CoverageUtilities.WEST);
    double e = regionMap.get(CoverageUtilities.EAST);
    double meanX = w + (e - w) / 2.0;
    double meanY = s + (n - s) / 2.0;
    Coordinate tmp = new Coordinate(meanX, meanY);
    MathTransform mathTransform = CRS.findMathTransform(inAspect.getCoordinateReferenceSystem(), DefaultGeographicCRS.WGS84);
    Coordinate newC = JTS.transform(tmp, null, mathTransform);
    avgLatitude = newC.y;
    RenderedImage idbasinImage = inBasins.getRenderedImage();
    idbasinImageIterator = RandomIterFactory.create(idbasinImage, null);
    RenderedImage elevImage = inElev.getRenderedImage();
    elevImageIterator = RandomIterFactory.create(elevImage, null);
    RenderedImage tmpImage = inCurvatures.getRenderedImage();
    curvatureImage = CoverageUtilities.createWritableRaster(tmpImage.getWidth(), tmpImage.getHeight(), null, null, null);
    RandomIter tmpIterator = RandomIterFactory.create(tmpImage, null);
    // TODO check what this is for?!?!?
    for (int i = 0; i < tmpImage.getHeight(); i++) {
        for (int j = 0; j < tmpImage.getWidth(); j++) {
            double value = tmpIterator.getSampleDouble(j, i, 0);
            curvatureImage.setSample(j, i, 0, value);
        }
    }
    RenderedImage aspectImage = inAspect.getRenderedImage();
    aspectImageIterator = RandomIterFactory.create(aspectImage, null);
    RenderedImage slopeImage = inSlope.getRenderedImage();
    slopeImageIterator = RandomIterFactory.create(slopeImage, null);
    avgLatitude *= (PI / 180.0);
    // $NON-NLS-1$
    pm.message(msg.message("eicalculator.preparing_inputs"));
    eibasinNum = prepareInputsOutputs();
    // $NON-NLS-1$
    pm.beginTask(msg.message("eicalculator.computing"), 6);
    for (int m = 0; m < 6; m++) {
        pm.worked(1);
        compute_EI(m + 1);
    }
    pm.done();
    average_EI(10, 6);
    // $NON-NLS-1$
    pm.beginTask(msg.message("eicalculator.calc_areas"), eibasinNum);
    for (int i = 0; i < eibasinNum; i++) {
        pm.worked(1);
        area(i);
    }
    pm.done();
    /*
         * putting the results together
         */
    outAltimetry = new ArrayList<EIAltimetry>();
    outEnergy = new ArrayList<EIEnergy>();
    outArea = new ArrayList<EIAreas>();
    for (int i = 0; i < eibasinNum; i++) {
        int realBasinId = index2idMap.get(i + 1);
        /*
             * ENERGY BANDS
             * 
             * Cycle over the virtual months:
             * 0: 22 DICEMBRE - 20 GENNAIO
             * 1: 21 GENNAIO - 20 FEBBRAIO
             * 2: 21 FEBBRAIO - 22 MARZO
             * 3: 23 MARZO - 22 APRILE
             * 4: 23 APRILE - 22 MAGGIO
             * 5: 23 MAGGIO - 22 GIUGNO
             */
        for (int j = 0; j < 6; j++) {
            for (int k = 0; k < pEi; k++) {
                EIEnergy tmpEi = new EIEnergy();
                // the basin id
                tmpEi.basinId = realBasinId;
                tmpEi.energeticBandId = k;
                tmpEi.virtualMonth = j;
                tmpEi.energyValue = eibasinEI[0][k][i];
                outEnergy.add(tmpEi);
            }
        }
        /*
             * ALTIMETRIC BANDS
             */
        for (int k = 0; k < pEs; k++) {
            EIAltimetry tmpAl = new EIAltimetry();
            tmpAl.basinId = realBasinId;
            tmpAl.altimetricBandId = k;
            tmpAl.elevationValue = eibasinES[k][i];
            tmpAl.bandRange = eibasinESrange[k][i];
            outAltimetry.add(tmpAl);
        }
        /*
             * AREAS
             */
        for (int j = 0; j < pEs; j++) {
            for (int k = 0; k < pEi; k++) {
                EIAreas tmpAr = new EIAreas();
                tmpAr.basinId = realBasinId;
                tmpAr.altimetricBandId = j;
                tmpAr.energyBandId = k;
                tmpAr.areaValue = eibasinA[j][k][i];
                outArea.add(tmpAr);
            }
        }
    }
}

Example 3 with EIAltimetry

use of org.hortonmachine.gears.io.eicalculator.EIAltimetry in project hortonmachine by TheHortonMachine.

the class TestEnergyIndexCalculator method testEnergyIndexCalculator.

/**
 * TODO make this test a bit more serious.
 *
 * @throws Exception
 */
public void testEnergyIndexCalculator() throws Exception {
    HashMap<String, Double> envelopeParams = HMTestMaps.getEnvelopeparams();
    CoordinateReferenceSystem crs = HMTestMaps.getCrs();
    // PrintStreamProgressMonitor pm = new PrintStreamProgressMonitor(System.out, System.out);
    double[][] aspectData = HMTestMaps.aspectDataRadiants;
    GridCoverage2D aspectCoverage = CoverageUtilities.buildCoverage("aspect", aspectData, envelopeParams, crs, true);
    double[][] nablaData = HMTestMaps.nablaData0;
    GridCoverage2D nablaCoverage = CoverageUtilities.buildCoverage("nabla", nablaData, envelopeParams, crs, true);
    double[][] pitData = HMTestMaps.pitData;
    GridCoverage2D pitCoverage = CoverageUtilities.buildCoverage("pit", pitData, envelopeParams, crs, true);
    double[][] slopeData = HMTestMaps.slopeData;
    GridCoverage2D slopeCoverage = CoverageUtilities.buildCoverage("slope", slopeData, envelopeParams, crs, true);
    double[][] subbasinsData = HMTestMaps.basinDataNN0;
    GridCoverage2D subbasinsCoverage = CoverageUtilities.buildCoverage("subbasins", subbasinsData, envelopeParams, crs, true);
    OmsEnergyIndexCalculator eiCalculator = new OmsEnergyIndexCalculator();
    eiCalculator.inAspect = aspectCoverage;
    eiCalculator.inCurvatures = nablaCoverage;
    eiCalculator.inElev = pitCoverage;
    eiCalculator.inSlope = slopeCoverage;
    eiCalculator.inBasins = subbasinsCoverage;
    eiCalculator.pDt = 1;
    eiCalculator.pEi = 2;
    eiCalculator.pEs = 2;
    eiCalculator.pm = new DummyProgressMonitor();
    eiCalculator.process();
    List<EIAltimetry> altimetricValues = eiCalculator.outAltimetry;
    List<EIEnergy> energeticValues = eiCalculator.outEnergy;
    List<EIAreas> areaValues = eiCalculator.outArea;
    EIAltimetry eiAltimetry = altimetricValues.get(0);
    assertEquals(1, eiAltimetry.basinId);
    assertEquals(0, eiAltimetry.altimetricBandId);
    assertEquals(737.5, eiAltimetry.elevationValue);
    assertEquals(75.0, eiAltimetry.bandRange);
    EIEnergy eiEnergy = energeticValues.get(0);
    assertEquals(1, eiEnergy.basinId);
    assertEquals(0, eiEnergy.energeticBandId);
    assertEquals(0, eiEnergy.virtualMonth);
    assertEquals(0.09808943859674346, eiEnergy.energyValue, 0.0001);
    EIAreas eiAreas = areaValues.get(0);
    assertEquals(1, eiAreas.basinId);
    assertEquals(0, eiAreas.altimetricBandId);
    assertEquals(0, eiAreas.energyBandId);
    assertEquals(0.0, eiAreas.areaValue, 0.0001);
}

Example 4 with EIAltimetry

use of org.hortonmachine.gears.io.eicalculator.EIAltimetry in project hortonmachine by TheHortonMachine.

the class EnergyIndexCalculator method process.

@Execute
public void process() throws Exception {
    OmsEnergyIndexCalculator energyindexcalculator = new OmsEnergyIndexCalculator();
    energyindexcalculator.inElev = getRaster(inElev);
    energyindexcalculator.inBasins = getRaster(inBasins);
    energyindexcalculator.inCurvatures = getRaster(inCurvatures);
    energyindexcalculator.inAspect = getRaster(inAspect);
    energyindexcalculator.inSlope = getRaster(inSlope);
    energyindexcalculator.pEs = pEs;
    energyindexcalculator.pEi = pEi;
    energyindexcalculator.pDt = pDt;
    energyindexcalculator.pm = pm;
    energyindexcalculator.process();
    List<EIAltimetry> outAltimetryObj = energyindexcalculator.outAltimetry;
    List<EIEnergy> outEnergyObj = energyindexcalculator.outEnergy;
    List<EIAreas> outAreaObj = energyindexcalculator.outArea;
    OmsEIAltimetryWriter altimetryWriter = new OmsEIAltimetryWriter();
    altimetryWriter.inAltimetry = outAltimetryObj;
    altimetryWriter.file = outAltimetry;
    altimetryWriter.write();
    altimetryWriter.close();
    OmsEIEnergyWriter energyWriter = new OmsEIEnergyWriter();
    energyWriter.inEnergy = outEnergyObj;
    energyWriter.file = outEnergy;
    energyWriter.write();
    energyWriter.close();
    OmsEIAreasWriter areasWriter = new OmsEIAreasWriter();
    areasWriter.inAreas = outAreaObj;
    areasWriter.file = outArea;
    areasWriter.write();
    areasWriter.close();
}

Example 5 with EIAltimetry

use of org.hortonmachine.gears.io.eicalculator.EIAltimetry in project hortonmachine by TheHortonMachine.

the class Jami method process.

@Execute
public void process() throws Exception {
    OmsEIAltimetryReader altim = new OmsEIAltimetryReader();
    altim.file = inAltimetry;
    altim.pSeparator = ",";
    altim.pm = pm;
    altim.read();
    List<EIAltimetry> altimList = altim.outAltimetry;
    altim.close();
    OmsEIAreasReader areas = new OmsEIAreasReader();
    areas.file = inAreas;
    areas.pSeparator = ",";
    areas.pm = pm;
    areas.read();
    List<EIAreas> areasList = areas.outAreas;
    areas.close();
    OmsTimeSeriesIteratorReader dataReader = new OmsTimeSeriesIteratorReader();
    dataReader.file = inMeteo;
    dataReader.fileNovalue = "-9999";
    dataReader.idfield = fStationid;
    // dataReader.tStart = "2005-05-01 00:00";
    // dataReader.tTimestep = 60;
    // dataReader.tEnd = "2005-05-01 03:00";
    dataReader.initProcess();
    OmsJami omsjami = new OmsJami();
    omsjami.inStations = getVector(inStations);
    omsjami.fStationid = fStationid;
    omsjami.fStationelev = fStationelev;
    omsjami.pBins = pBins;
    omsjami.pNum = pNum;
    omsjami.inInterpolate = getVector(inInterpolate);
    omsjami.fBasinid = fBasinid;
    omsjami.pType = pType;
    omsjami.defaultRh = defaultRh;
    omsjami.defaultW = defaultW;
    omsjami.pHtmin = pHtmin;
    omsjami.pHtmax = pHtmax;
    omsjami.defaultDtday = defaultDtday;
    omsjami.defaultDtmonth = defaultDtmonth;
    omsjami.defaultTolltmin = defaultTolltmin;
    omsjami.defaultTolltmax = defaultTolltmax;
    omsjami.inAltimetry = altimList;
    omsjami.inAreas = areasList;
    omsjami.pm = pm;
    omsjami.doProcess = doProcess;
    omsjami.doReset = doReset;
    OmsTimeSeriesIteratorWriter writerInterpolated = null;
    OmsTimeSeriesIteratorWriter writerInterpolatedBand = null;
    pm.beginTask("Processing...", IHMProgressMonitor.UNKNOWN);
    while (dataReader.doProcess) {
        dataReader.nextRecord();
        if (writerInterpolated == null && outInterpolated != null) {
            writerInterpolated = new OmsTimeSeriesIteratorWriter();
            writerInterpolated.file = outInterpolated;
            writerInterpolated.tStart = dataReader.tStart;
            writerInterpolated.tTimestep = dataReader.tTimestep;
        }
        if (writerInterpolatedBand == null && outInterpolatedBand != null) {
            writerInterpolatedBand = new OmsTimeSeriesIteratorWriter();
            writerInterpolatedBand.file = outInterpolatedBand;
            writerInterpolatedBand.tStart = dataReader.tStart;
            writerInterpolatedBand.tTimestep = dataReader.tTimestep;
        }
        pm.message("timestep: " + dataReader.tCurrent);
        omsjami.tCurrent = dataReader.tCurrent;
        HashMap<Integer, double[]> id2ValueMap = dataReader.outData;
        omsjami.inMeteo = id2ValueMap;
        omsjami.process();
        // write csv data
        writerInterpolated.inData = omsjami.outInterpolated;
        writerInterpolated.writeNextLine();
        writerInterpolatedBand.inData = omsjami.outInterpolatedBand;
        writerInterpolatedBand.writeNextLine();
    }
    pm.done();
    dataReader.close();
    if (writerInterpolated != null)
        writerInterpolated.close();
    if (writerInterpolatedBand != null)
        writerInterpolatedBand.close();
}