Examples with RandomIter javax.media.jai.iterator.RandomIter used on opensource projects

Example 66 with RandomIter

use of javax.media.jai.iterator.RandomIter in project hortonmachine by TheHortonMachine.

the class OmsHillshade method normalVector.

protected WritableRaster normalVector(WritableRaster pitWR, double res) {
    int minX = pitWR.getMinX();
    int minY = pitWR.getMinY();
    int rows = pitWR.getHeight();
    int cols = pitWR.getWidth();
    RandomIter pitIter = RandomIterFactory.create(pitWR, null);
    /*
         * Initialize the Image of the normal vector in the central point of the
         * cells, which have 3 components so the Image have 3 bands..
         */
    SampleModel sm = RasterFactory.createBandedSampleModel(5, cols, rows, 3);
    WritableRaster tmpNormalVectorWR = CoverageUtilities.createWritableRaster(cols, rows, null, sm, 0.0);
    WritableRandomIter tmpNormaIter = RandomIterFactory.createWritable(tmpNormalVectorWR, null);
    /*
         * apply the corripio's formula (is the formula (3) in the article)
         */
    for (int j = minY; j < minX + rows - 1; j++) {
        for (int i = minX; i < minX + cols - 1; i++) {
            double zij = pitIter.getSampleDouble(i, j, 0);
            double zidxj = pitIter.getSampleDouble(i + 1, j, 0);
            double zijdy = pitIter.getSampleDouble(i, j + 1, 0);
            double zidxjdy = pitIter.getSampleDouble(i + 1, j + 1, 0);
            double firstComponent = 0.5 * res * (zij - zidxj + zijdy - zidxjdy);
            double secondComponent = 0.5 * res * (zij + zidxj - zijdy - zidxjdy);
            double thirthComponent = (res * res);
            double den = Math.sqrt(firstComponent * firstComponent + secondComponent * secondComponent + thirthComponent * thirthComponent);
            tmpNormaIter.setPixel(i, j, new double[] { firstComponent / den, secondComponent / den, thirthComponent / den });
        }
    }
    pitIter.done();
    return tmpNormalVectorWR;
}

Example 67 with RandomIter

use of javax.media.jai.iterator.RandomIter in project hortonmachine by TheHortonMachine.

the class OmsDebrisFlow method process.

@Execute
public void process() throws Exception {
    RegionMap regionMap = CoverageUtilities.getRegionParamsFromGridCoverage(inElev);
    int cols = regionMap.getCols();
    int rows = regionMap.getRows();
    double xRes = regionMap.getXres();
    double yRes = regionMap.getYres();
    double west = regionMap.getWest();
    double east = regionMap.getEast();
    double south = regionMap.getSouth();
    double north = regionMap.getNorth();
    if (!isBetween(pEasting, west, east) || !isBetween(pNorthing, south, north)) {
        throw new ModelsIllegalargumentException("Input coordinates have to be within the map boundaries.", this, pm);
    }
    double thresArea = pMcoeff * pow(pVolume, (2.0 / 3.0));
    GridGeometry2D gridGeometry = inElev.getGridGeometry();
    int[] colRow = CoverageUtilities.colRowFromCoordinate(new Coordinate(pEasting, pNorthing), gridGeometry, null);
    RandomIter elevIter = CoverageUtilities.getRandomIterator(inElev);
    int startCol = colRow[0];
    int startRow = colRow[1];
    double startValue = elevIter.getSampleDouble(startCol, startRow, 0);
    if (isNovalue(startValue)) {
        throw new ModelsIllegalargumentException("Input coordinates are on a novalue elevation point.", this, pm);
    }
    WritableRaster mcsWR = CoverageUtilities.createWritableRaster(cols, rows, null, null, HMConstants.doubleNovalue);
    WritableRandomIter probIter = RandomIterFactory.createWritable(mcsWR, null);
    Random flatRnd = new Random();
    int processedMc = 0;
    for (int mc = 0; mc < pMontecarlo; mc++) {
        pm.message("Montecarlo n." + mc);
        processedMc = mc;
        /*
             * for every Montecarlo loop, get a flow path
             */
        int centerCol = startCol;
        int centerRow = startRow;
        TreeSet<Point> touchedPoints = new TreeSet<Point>();
        touchedPoints.add(new Point(centerCol, centerRow));
        boolean doStop = false;
        // cicle over every cell neighbours along the way
        // -System.currentTimeMillis());
        Random randomGenerator = new Random();
        do {
            // System.out.println(centerCol + "/" + centerRow + " --- " + cols + "/" + rows);
            double centerValue = elevIter.getSampleDouble(centerCol, centerRow, 0);
            List<SlopeProbability> spList = new ArrayList<SlopeProbability>();
            double slopeSum = 0;
            for (int x = -1; x <= 1; x++) {
                for (int y = -1; y <= 1; y++) {
                    if (x == 0 && y == 0) {
                        continue;
                    }
                    int tmpCol = centerCol + x;
                    int tmpRow = centerRow + y;
                    if (touchedPoints.contains(new Point(tmpCol, tmpRow))) {
                        continue;
                    }
                    // if point is outside jump it
                    if (!isBetween(tmpCol, 0, cols - 1) || !isBetween(tmpRow, 0, rows - 1)) {
                        continue;
                    }
                    // if point is novalue, jump it
                    double nextValue = elevIter.getSampleDouble(tmpCol, tmpRow, 0);
                    if (isNovalue(nextValue)) {
                        continue;
                    }
                    double distance = pythagoras(abs((tmpCol - centerCol) * xRes), abs((tmpRow - centerRow) * yRes));
                    double slope = (nextValue - centerValue) / distance;
                    // we take only negative and 0 slope, downhill
                    if (slope > 0) {
                        continue;
                    }
                    slope = abs(slope);
                    SlopeProbability sp = new SlopeProbability();
                    sp.fromCol = centerCol;
                    sp.fromRow = centerRow;
                    sp.fromElev = centerValue;
                    sp.toCol = tmpCol;
                    sp.toRow = tmpRow;
                    sp.toElev = nextValue;
                    sp.slope = slope;
                    slopeSum = slopeSum + slope;
                    spList.add(sp);
                }
            }
            if (spList.size() == 0) {
                /*
                     * touched border or slope is not negative
                     */
                doStop = true;
            } else {
                // get a random number between 0 and 1
                double random = randomGenerator.nextDouble();
                if (spList.size() == 1) {
                    // direction is only one
                    SlopeProbability sp = spList.get(0);
                    centerCol = sp.toCol;
                    centerRow = sp.toRow;
                } else {
                    Collections.sort(spList);
                    /*
                         * case in which the slopes are all 0
                         */
                    if (dEq(slopeSum, 0.0)) {
                        // choose a random and go on
                        int size = spList.size();
                        double rnd = flatRnd.nextDouble();
                        int index = (int) round(rnd * size) - 1;
                        if (index < 0)
                            index = 0;
                        SlopeProbability sp = spList.get(index);
                        centerCol = sp.toCol;
                        centerRow = sp.toRow;
                    } else /*
                         * normal case in which the slopes have a value
                         */
                    {
                        // cumulate the probability
                        for (int i = 0; i < spList.size(); i++) {
                            SlopeProbability sp = spList.get(i);
                            double p = sp.slope / slopeSum;
                            sp.probability = p;
                            if (i != 0) {
                                SlopeProbability tmpSp = spList.get(i - 1);
                                sp.probability = sp.probability + tmpSp.probability;
                            }
                        }
                        for (int i = 1; i < spList.size(); i++) {
                            SlopeProbability sp1 = spList.get(i - 1);
                            SlopeProbability sp2 = spList.get(i);
                            // if (random < sp1.probability) {
                            if (random < sp1.probability) {
                                centerCol = sp1.toCol;
                                centerRow = sp1.toRow;
                                break;
                            // } else if (random >= sp1.probability && random <
                            // sp2.probability)
                            // {
                            } else if (random >= sp1.probability && random < sp2.probability) {
                                centerCol = sp2.toCol;
                                centerRow = sp2.toRow;
                                break;
                            }
                        }
                    }
                }
                touchedPoints.add(new Point(centerCol, centerRow));
                double outValue = probIter.getSampleDouble(centerCol, centerRow, 0);
                if (isNovalue(outValue)) {
                    outValue = 0.0;
                }
                probIter.setSample(centerCol, centerRow, 0, outValue + 1.0);
            }
        } while (!doStop);
        /*
             * check if the max area is flooded
             */
        int floodedCellNum = 0;
        for (int r = 0; r < rows; r++) {
            for (int c = 0; c < cols; c++) {
                double value = probIter.getSampleDouble(c, r, 0);
                if (isNovalue(value)) {
                    continue;
                }
                floodedCellNum++;
            }
        }
        double floodedArea = floodedCellNum * xRes * yRes;
        if (thresArea <= floodedArea) {
            break;
        }
    }
    double probSum = 0.0;
    double validCells = 0.0;
    for (int r = 0; r < rows; r++) {
        for (int c = 0; c < cols; c++) {
            double prob = probIter.getSampleDouble(c, r, 0);
            if (isNovalue(prob)) {
                continue;
            }
            double newProb = prob / (processedMc - 1);
            probIter.setSample(c, r, 0, newProb);
            probSum = probSum + sqrt(newProb);
            validCells++;
        }
    }
    /*
         * calculate deposition
         */
    double avgProb = probSum / validCells;
    double avgHeight = pDcoeff * pow(pVolume, 1.0 / 3.0);
    WritableRaster depoWR = CoverageUtilities.createWritableRaster(cols, rows, null, null, HMConstants.doubleNovalue);
    WritableRandomIter depoIter = RandomIterFactory.createWritable(depoWR, null);
    for (int r = 0; r < rows; r++) {
        for (int c = 0; c < cols; c++) {
            double probValue = probIter.getSampleDouble(c, r, 0);
            if (isNovalue(probValue)) {
                continue;
            }
            double depoValue = avgHeight * sqrt(probValue) / avgProb;
            depoIter.setSample(c, r, 0, depoValue);
        }
    }
    outMcs = CoverageUtilities.buildCoverage("mcs", mcsWR, regionMap, inElev.getCoordinateReferenceSystem());
    outDepo = CoverageUtilities.buildCoverage("depo", depoWR, regionMap, inElev.getCoordinateReferenceSystem());
}

Example 68 with RandomIter

use of javax.media.jai.iterator.RandomIter in project hortonmachine by TheHortonMachine.

the class OmsDebrisTriggerCnr method process.

@Execute
public void process() throws Exception {
    checkNull(inElev, inNet, inTca);
    // calculate gradient map degrees
    OmsGradient gradient = new OmsGradient();
    gradient.inElev = inElev;
    gradient.pMode = Variables.FINITE_DIFFERENCES;
    gradient.doDegrees = true;
    gradient.pm = pm;
    gradient.process();
    GridCoverage2D gradientCoverageDeg = gradient.outSlope;
    // calculate gradient map %
    gradient = new OmsGradient();
    gradient.inElev = inElev;
    gradient.pMode = Variables.FINITE_DIFFERENCES;
    gradient.doDegrees = false;
    gradient.pm = pm;
    gradient.process();
    GridCoverage2D gradientCoverageTan = gradient.outSlope;
    // ritaglio della mappa di gradient lungo il reticolo
    // idrografico ed estrazione delle sole celle con
    // * pendenza minore di 38 gradi
    // * area cumulata minore di 10 km2
    RegionMap regionMap = CoverageUtilities.getRegionParamsFromGridCoverage(inElev);
    int cols = regionMap.getCols();
    int rows = regionMap.getRows();
    double xres = regionMap.getXres();
    double yres = regionMap.getYres();
    RenderedImage netRI = inNet.getRenderedImage();
    RandomIter netIter = RandomIterFactory.create(netRI, null);
    RenderedImage tcaRI = inTca.getRenderedImage();
    RandomIter tcaIter = RandomIterFactory.create(tcaRI, null);
    RenderedImage gradientDegRI = gradientCoverageDeg.getRenderedImage();
    RandomIter gradientDegIter = RandomIterFactory.create(gradientDegRI, null);
    RenderedImage gradientTanRI = gradientCoverageTan.getRenderedImage();
    RandomIter gradientTanIter = RandomIterFactory.create(gradientTanRI, null);
    WritableRaster outputWR = CoverageUtilities.createWritableRaster(cols, rows, null, null, HMConstants.doubleNovalue);
    WritableRandomIter outputIter = RandomIterFactory.createWritable(outputWR, null);
    pm.beginTask("Extracting trigger points...", cols);
    for (int r = 0; r < rows; r++) {
        for (int c = 0; c < cols; c++) {
            double net = netIter.getSampleDouble(c, r, 0);
            // all only along the network
            if (!isNovalue(net)) {
                double tca = tcaIter.getSampleDouble(c, r, 0);
                // tca in km2 along the net
                double tcaKm2 = tca * xres * yres / 1000000;
                // gradient in degrees along the net
                double gradientDeg = gradientDegIter.getSampleDouble(c, r, 0);
                // gradient in tan along the net
                double gradientTan = gradientTanIter.getSampleDouble(c, r, 0);
                /*
                     * calculate the trigger threshold:
                     * 
                     *  S = 0.32 * A^-0.2
                     *  where:
                     *   S = gradient in m/m
                     *   A = tca in km2
                     */
                double triggerThreshold = 0.32 * pow(tcaKm2, -0.2);
                if (// 
                gradientTan > triggerThreshold && // 
                gradientDeg < pGradthres && tcaKm2 < pTcathres) {
                    // we have a trigger point
                    outputIter.setSample(c, r, 0, triggerThreshold);
                }
            }
        }
        pm.worked(1);
    }
    pm.done();
    outTriggers = CoverageUtilities.buildCoverage("triggers", outputWR, regionMap, inElev.getCoordinateReferenceSystem());
}

Example 69 with RandomIter

use of javax.media.jai.iterator.RandomIter in project hortonmachine by TheHortonMachine.

the class OmsInsolation method normalVector.

protected WritableRaster normalVector(WritableRaster pitWR, double res) {
    int minX = pitWR.getMinX();
    int minY = pitWR.getMinY();
    int rows = pitWR.getHeight();
    int cols = pitWR.getWidth();
    RandomIter pitIter = RandomIterFactory.create(pitWR, null);
    /*
         * Initializa the Image of the normal vector in the central point of the
         * cells, which have 3 components so the Image have 3 bands..
         */
    SampleModel sm = RasterFactory.createBandedSampleModel(5, cols, rows, 3);
    WritableRaster tmpNormalVectorWR = CoverageUtilities.createWritableRaster(cols, rows, null, sm, 0.0);
    WritableRandomIter tmpNormalIter = RandomIterFactory.createWritable(tmpNormalVectorWR, null);
    /*
         * apply the corripio's formula (is the formula (3) in the article)
         */
    for (int j = minY; j < minX + rows - 1; j++) {
        for (int i = minX; i < minX + cols - 1; i++) {
            double zij = pitIter.getSampleDouble(i, j, 0);
            double zidxj = pitIter.getSampleDouble(i + 1, j, 0);
            double zijdy = pitIter.getSampleDouble(i, j + 1, 0);
            double zidxjdy = pitIter.getSampleDouble(i + 1, j + 1, 0);
            double firstComponent = res * (zij - zidxj + zijdy - zidxjdy);
            double secondComponent = res * (zij + zidxj - zijdy - zidxjdy);
            double thirthComponent = 2 * (res * res);
            double den = Math.sqrt(firstComponent * firstComponent + secondComponent * secondComponent + thirthComponent * thirthComponent);
            tmpNormalIter.setPixel(i, j, new double[] { firstComponent / den, secondComponent / den, thirthComponent / den });
        }
    }
    pitIter.done();
    return tmpNormalVectorWR;
}

Example 70 with RandomIter

use of javax.media.jai.iterator.RandomIter in project hortonmachine by TheHortonMachine.

the class OmsHazardClassifier method process.

@Execute
public void process() throws Exception {
    checkNull(inIntensityTr100, inIntensityTr200, inIntensityTr30);
    RegionMap regionMap = CoverageUtilities.getRegionParamsFromGridCoverage(inIntensityTr30);
    int nCols = regionMap.getCols();
    int nRows = regionMap.getRows();
    RandomIter tr200Iter = CoverageUtilities.getRandomIterator(inIntensityTr200);
    RandomIter tr100Iter = CoverageUtilities.getRandomIterator(inIntensityTr100);
    RandomIter tr30Iter = CoverageUtilities.getRandomIterator(inIntensityTr30);
    WritableRaster outIP1WR = CoverageUtilities.createWritableRaster(nCols, nRows, null, null, doubleNovalue);
    WritableRandomIter outIP1Iter = RandomIterFactory.createWritable(outIP1WR, null);
    WritableRaster outIP2WR = CoverageUtilities.createWritableRaster(nCols, nRows, null, null, doubleNovalue);
    WritableRandomIter outIP2Iter = RandomIterFactory.createWritable(outIP2WR, null);
    pm.beginTask("Processing map...", nRows);
    for (int r = 0; r < nRows; r++) {
        if (isCanceled(pm)) {
            return;
        }
        for (int c = 0; c < nCols; c++) {
            double tr30 = tr30Iter.getSampleDouble(c, r, 0);
            double tr100 = tr100Iter.getSampleDouble(c, r, 0);
            double tr200 = tr200Iter.getSampleDouble(c, r, 0);
            if (isNovalue(tr30) && isNovalue(tr100) && isNovalue(tr200)) {
                continue;
            }
            double tmpTr30;
            if (isNovalue(tr30)) {
                tmpTr30 = Double.NEGATIVE_INFINITY;
            } else if (dEq(tr30, INTENSITY_LOW)) {
                tmpTr30 = 3.0;
            } else if (dEq(tr30, INTENSITY_MEDIUM)) {
                tmpTr30 = 6.0;
            } else if (dEq(tr30, INTENSITY_HIGH)) {
                tmpTr30 = 9.0;
            } else {
                throw new ModelsIllegalargumentException("Unknown tr30 value: " + tr30, this, pm);
            }
            double tmpTr100;
            if (isNovalue(tr100)) {
                tmpTr100 = Double.NEGATIVE_INFINITY;
            } else if (dEq(tr100, INTENSITY_LOW)) {
                tmpTr100 = 2.0;
            } else if (dEq(tr100, INTENSITY_MEDIUM)) {
                tmpTr100 = 5.0;
            } else if (dEq(tr100, INTENSITY_HIGH)) {
                tmpTr100 = 8.0;
            } else {
                throw new ModelsIllegalargumentException("Unknown tr100 value: " + tr100, this, pm);
            }
            double tmpTr200;
            if (isNovalue(tr200)) {
                tmpTr200 = Double.NEGATIVE_INFINITY;
            } else if (dEq(tr200, INTENSITY_LOW)) {
                tmpTr200 = 1.0;
            } else if (dEq(tr200, INTENSITY_MEDIUM)) {
                tmpTr200 = 4.0;
            } else if (dEq(tr200, INTENSITY_HIGH)) {
                tmpTr200 = 7.0;
            } else {
                throw new ModelsIllegalargumentException("Unknown tr200 value: " + tr200, this, pm);
            }
            int maxValue = (int) max(tmpTr30, max(tmpTr100, tmpTr200));
            double[] reclassIP1 = { // 
            Double.NaN, // 1
            2, // 2
            3, // 3
            3, // 4
            3, // 5
            3, // 6
            4, // 7
            4, // 8
            4, // 9
            4 };
            double[] reclassIP2 = { // 
            Double.NaN, // 1
            2, // 2
            2, // 3
            3, // 4
            3, // 5
            3, // 6
            3, // 7
            4, // 8
            4, // 9
            4 };
            if (maxValue < 1 || maxValue > (reclassIP1.length - 1)) {
                throw new ModelsIllegalargumentException("Unknown max value from tr30/100/200: " + maxValue, this, pm);
            }
            double ip1 = reclassIP1[(int) maxValue];
            double ip2 = reclassIP2[(int) maxValue];
            outIP1Iter.setSample(c, r, 0, ip1);
            outIP2Iter.setSample(c, r, 0, ip2);
        }
        pm.worked(1);
    }
    pm.done();
    outHazardIP1 = CoverageUtilities.buildCoverage("ip1", outIP1WR, regionMap, inIntensityTr100.getCoordinateReferenceSystem());
    outHazardIP2 = CoverageUtilities.buildCoverage("ip2", outIP2WR, regionMap, inIntensityTr100.getCoordinateReferenceSystem());
}