Search in sources :

Example 1 with DemTile

use of org.jlinda.core.geom.DemTile in project s1tbx by senbox-org.

the class InterferogramOp method computeTileStackForNormalProduct.

private void computeTileStackForNormalProduct(final Map<Band, Tile> targetTileMap, Rectangle targetRectangle, final ProgressMonitor pm) throws OperatorException {
    try {
        final BorderExtender border = BorderExtender.createInstance(BorderExtender.BORDER_ZERO);
        final int y0 = targetRectangle.y;
        final int yN = y0 + targetRectangle.height - 1;
        final int x0 = targetRectangle.x;
        final int xN = targetRectangle.x + targetRectangle.width - 1;
        final Window tileWindow = new Window(y0, yN, x0, xN);
        DemTile demTile = null;
        if (subtractTopographicPhase) {
            demTile = TopoPhase.getDEMTile(tileWindow, targetMap, dem, demNoDataValue, demSamplingLat, demSamplingLon, tileExtensionPercent);
            if (demTile.getData().length < 3 || demTile.getData()[0].length < 3) {
                throw new OperatorException("The resolution of the selected DEM is too low, " + "please select DEM with higher resolution.");
            }
        }
        // parameters for coherence calculation
        final int cohx0 = targetRectangle.x - (cohWinRg - 1) / 2;
        final int cohy0 = targetRectangle.y - (cohWinAz - 1) / 2;
        final int cohw = targetRectangle.width + cohWinRg - 1;
        final int cohh = targetRectangle.height + cohWinAz - 1;
        final Rectangle rect = new Rectangle(cohx0, cohy0, cohw, cohh);
        final Window cohTileWindow = new Window(cohy0, cohy0 + cohh - 1, cohx0, cohx0 + cohw - 1);
        DemTile cohDemTile = null;
        if (subtractTopographicPhase) {
            cohDemTile = TopoPhase.getDEMTile(cohTileWindow, targetMap, dem, demNoDataValue, demSamplingLat, demSamplingLon, tileExtensionPercent);
        }
        for (String ifgKey : targetMap.keySet()) {
            final ProductContainer product = targetMap.get(ifgKey);
            final Tile mstTileReal = getSourceTile(product.sourceMaster.realBand, targetRectangle, border);
            final Tile mstTileImag = getSourceTile(product.sourceMaster.imagBand, targetRectangle, border);
            final ComplexDoubleMatrix dataMaster = TileUtilsDoris.pullComplexDoubleMatrix(mstTileReal, mstTileImag);
            final Tile slvTileReal = getSourceTile(product.sourceSlave.realBand, targetRectangle, border);
            final Tile slvTileImag = getSourceTile(product.sourceSlave.imagBand, targetRectangle, border);
            final ComplexDoubleMatrix dataSlave = TileUtilsDoris.pullComplexDoubleMatrix(slvTileReal, slvTileImag);
            if (subtractFlatEarthPhase) {
                final DoubleMatrix flatEarthPhase = computeFlatEarthPhase(x0, xN, dataMaster.columns, y0, yN, dataMaster.rows, 0, sourceImageWidth - 1, 0, sourceImageHeight - 1, product.sourceSlave.name);
                final ComplexDoubleMatrix complexReferencePhase = new ComplexDoubleMatrix(MatrixFunctions.cos(flatEarthPhase), MatrixFunctions.sin(flatEarthPhase));
                dataSlave.muli(complexReferencePhase);
                if (OUTPUT_PHASE) {
                    saveFlatEarthPhase(x0, xN, y0, yN, flatEarthPhase, product, targetTileMap);
                }
            }
            if (subtractTopographicPhase) {
                final TopoPhase topoPhase = TopoPhase.computeTopoPhase(product, tileWindow, demTile, outputElevation, false);
                final ComplexDoubleMatrix ComplexTopoPhase = new ComplexDoubleMatrix(MatrixFunctions.cos(new DoubleMatrix(topoPhase.demPhase)), MatrixFunctions.sin(new DoubleMatrix(topoPhase.demPhase)));
                dataSlave.muli(ComplexTopoPhase);
                if (OUTPUT_PHASE) {
                    saveTopoPhase(x0, xN, y0, yN, topoPhase.demPhase, product, targetTileMap);
                }
                if (outputElevation) {
                    saveElevation(x0, xN, y0, yN, topoPhase.elevation, product, targetTileMap);
                }
                if (outputLatLon) {
                    final TopoPhase topoPhase1 = TopoPhase.computeTopoPhase(product, tileWindow, demTile, false, true);
                    saveLatLon(x0, xN, y0, yN, topoPhase1.latitude, topoPhase1.longitude, product, targetTileMap);
                }
            }
            dataMaster.muli(dataSlave.conji());
            saveInterferogram(dataMaster, product, targetTileMap, targetRectangle);
            // coherence calculation
            if (includeCoherence) {
                final Tile mstTileReal2 = getSourceTile(product.sourceMaster.realBand, rect, border);
                final Tile mstTileImag2 = getSourceTile(product.sourceMaster.imagBand, rect, border);
                final Tile slvTileReal2 = getSourceTile(product.sourceSlave.realBand, rect, border);
                final Tile slvTileImag2 = getSourceTile(product.sourceSlave.imagBand, rect, border);
                final ComplexDoubleMatrix dataMaster2 = TileUtilsDoris.pullComplexDoubleMatrix(mstTileReal2, mstTileImag2);
                final ComplexDoubleMatrix dataSlave2 = TileUtilsDoris.pullComplexDoubleMatrix(slvTileReal2, slvTileImag2);
                if (subtractFlatEarthPhase) {
                    final DoubleMatrix flatEarthPhase = computeFlatEarthPhase(cohx0, cohx0 + cohw - 1, cohw, cohy0, cohy0 + cohh - 1, cohh, 0, sourceImageWidth - 1, 0, sourceImageHeight - 1, product.sourceSlave.name);
                    final ComplexDoubleMatrix complexReferencePhase = new ComplexDoubleMatrix(MatrixFunctions.cos(flatEarthPhase), MatrixFunctions.sin(flatEarthPhase));
                    dataSlave2.muli(complexReferencePhase);
                }
                if (subtractTopographicPhase) {
                    final TopoPhase topoPhase = TopoPhase.computeTopoPhase(product, cohTileWindow, cohDemTile, false);
                    final ComplexDoubleMatrix ComplexTopoPhase = new ComplexDoubleMatrix(MatrixFunctions.cos(new DoubleMatrix(topoPhase.demPhase)), MatrixFunctions.sin(new DoubleMatrix(topoPhase.demPhase)));
                    dataSlave2.muli(ComplexTopoPhase);
                }
                for (int i = 0; i < dataMaster2.length; i++) {
                    double tmp = norm(dataMaster2.get(i));
                    dataMaster2.put(i, dataMaster2.get(i).mul(dataSlave2.get(i).conj()));
                    dataSlave2.put(i, new ComplexDouble(norm(dataSlave2.get(i)), tmp));
                }
                DoubleMatrix cohMatrix = SarUtils.coherence2(dataMaster2, dataSlave2, cohWinAz, cohWinRg);
                saveCoherence(cohMatrix, product, targetTileMap, targetRectangle);
            }
        }
    } catch (Throwable e) {
        OperatorUtils.catchOperatorException(getId(), e);
    } finally {
        pm.done();
    }
}
Also used : Window(org.jlinda.core.Window) BorderExtender(javax.media.jai.BorderExtender) DemTile(org.jlinda.core.geom.DemTile) Tile(org.esa.snap.core.gpf.Tile) Point(org.jlinda.core.Point) DemTile(org.jlinda.core.geom.DemTile) TopoPhase(org.jlinda.core.geom.TopoPhase) OperatorException(org.esa.snap.core.gpf.OperatorException)

Example 2 with DemTile

use of org.jlinda.core.geom.DemTile in project s1tbx by senbox-org.

the class CoherenceOp method computePartialTile.

private void computePartialTile(final int subSwathIndex, final int burstIndex, final int firstLineIdx, final Rectangle targetRectangle, final Map<Band, Tile> targetTileMap) throws OperatorException {
    try {
        final BorderExtender border = BorderExtender.createInstance(BorderExtender.BORDER_ZERO);
        final int y0 = targetRectangle.y;
        final int yN = y0 + targetRectangle.height - 1;
        final int x0 = targetRectangle.x;
        final int xN = x0 + targetRectangle.width - 1;
        final int cohx0 = targetRectangle.x - (cohWinRg - 1) / 2;
        final int cohy0 = targetRectangle.y - (cohWinAz - 1) / 2;
        final int cohw = targetRectangle.width + cohWinRg - 1;
        final int cohh = targetRectangle.height + cohWinAz - 1;
        final Rectangle extRect = new Rectangle(cohx0, cohy0, cohw, cohh);
        final org.jlinda.core.Window tileWindow = new org.jlinda.core.Window(cohy0 - firstLineIdx, cohy0 + cohh - 1 - firstLineIdx, cohx0, cohx0 + cohw - 1);
        final SLCImage mstMeta = targetMap.values().iterator().next().sourceMaster.metaData.clone();
        updateMstMetaData(burstIndex, mstMeta);
        final Orbit mstOrbit = targetMap.values().iterator().next().sourceMaster.orbit;
        DemTile demTile = null;
        if (subtractTopographicPhase) {
            demTile = TopoPhase.getDEMTile(tileWindow, mstMeta, mstOrbit, dem, demNoDataValue, demSamplingLat, demSamplingLon, tileExtensionPercent);
            if (demTile == null) {
                throw new OperatorException("The selected DEM has no overlap with the image or is invalid.");
            }
            if (demTile.getData().length < 3 || demTile.getData()[0].length < 3) {
                throw new OperatorException("The resolution of the selected DEM is too low, " + "please select DEM with higher resolution.");
            }
        }
        final int minLine = 0;
        final int maxLine = subSwath[subSwathIndex - 1].linesPerBurst - 1;
        final int minPixel = 0;
        final int maxPixel = subSwath[subSwathIndex - 1].samplesPerBurst - 1;
        for (String cohKey : targetMap.keySet()) {
            final ProductContainer product = targetMap.get(cohKey);
            final SLCImage slvMeta = product.sourceSlave.metaData.clone();
            updateSlvMetaData(product, burstIndex, slvMeta);
            final Orbit slvOrbit = product.sourceSlave.orbit;
            final Tile mstTileReal = getSourceTile(product.sourceMaster.realBand, extRect, border);
            final Tile mstTileImag = getSourceTile(product.sourceMaster.imagBand, extRect, border);
            final ComplexDoubleMatrix dataMaster = TileUtilsDoris.pullComplexDoubleMatrix(mstTileReal, mstTileImag);
            final Tile slvTileReal = getSourceTile(product.sourceSlave.realBand, extRect, border);
            final Tile slvTileImag = getSourceTile(product.sourceSlave.imagBand, extRect, border);
            final ComplexDoubleMatrix dataSlave = TileUtilsDoris.pullComplexDoubleMatrix(slvTileReal, slvTileImag);
            final String polynomialName = product.sourceSlave.name + '_' + (subSwathIndex - 1) + '_' + burstIndex;
            if (subtractFlatEarthPhase) {
                final DoubleMatrix flatEarthPhase = computeFlatEarthPhase(cohx0, cohx0 + cohw - 1, cohw, cohy0 - firstLineIdx, cohy0 + cohh - 1 - firstLineIdx, cohh, minPixel, maxPixel, minLine, maxLine, polynomialName);
                final ComplexDoubleMatrix complexReferencePhase = new ComplexDoubleMatrix(MatrixFunctions.cos(flatEarthPhase), MatrixFunctions.sin(flatEarthPhase));
                dataSlave.muli(complexReferencePhase);
                if (OUTPUT_PHASE) {
                    saveFlatEarthPhase(x0, xN, y0, yN, flatEarthPhase, product, targetTileMap);
                }
            }
            if (subtractTopographicPhase) {
                TopoPhase topoPhase = TopoPhase.computeTopoPhase(mstMeta, mstOrbit, slvMeta, slvOrbit, tileWindow, demTile, false);
                final ComplexDoubleMatrix ComplexTopoPhase = new ComplexDoubleMatrix(MatrixFunctions.cos(new DoubleMatrix(topoPhase.demPhase)), MatrixFunctions.sin(new DoubleMatrix(topoPhase.demPhase)));
                dataSlave.muli(ComplexTopoPhase);
                if (OUTPUT_PHASE) {
                    saveTopoPhase(x0, xN, y0, yN, topoPhase.demPhase, product, targetTileMap);
                }
            }
            for (int i = 0; i < dataMaster.length; i++) {
                double tmp = norm(dataMaster.get(i));
                dataMaster.put(i, dataMaster.get(i).mul(dataSlave.get(i).conj()));
                dataSlave.put(i, new ComplexDouble(norm(dataSlave.get(i)), tmp));
            }
            DoubleMatrix cohMatrix = SarUtils.coherence2(dataMaster, dataSlave, cohWinAz, cohWinRg);
            saveCoherence(cohMatrix, product, targetTileMap, targetRectangle);
        }
    } catch (Throwable e) {
        OperatorUtils.catchOperatorException(getId(), e);
    }
}
Also used : Orbit(org.jlinda.core.Orbit) BorderExtender(javax.media.jai.BorderExtender) DemTile(org.jlinda.core.geom.DemTile) Tile(org.esa.snap.core.gpf.Tile) Point(org.jlinda.core.Point) GeoPoint(org.jlinda.core.GeoPoint) DoubleMatrix(org.jblas.DoubleMatrix) ComplexDoubleMatrix(org.jblas.ComplexDoubleMatrix) ComplexDouble(org.jblas.ComplexDouble) DemTile(org.jlinda.core.geom.DemTile) TopoPhase(org.jlinda.core.geom.TopoPhase) SLCImage(org.jlinda.core.SLCImage) OperatorException(org.esa.snap.core.gpf.OperatorException) ComplexDoubleMatrix(org.jblas.ComplexDoubleMatrix)

Example 3 with DemTile

use of org.jlinda.core.geom.DemTile in project s1tbx by senbox-org.

the class CoherenceOp method computeTileForNormalProduct.

private void computeTileForNormalProduct(final Map<Band, Tile> targetTileMap, Rectangle targetRectangle, ProgressMonitor pm) throws OperatorException {
    try {
        final BorderExtender border = BorderExtender.createInstance(BorderExtender.BORDER_ZERO);
        final int y0 = targetRectangle.y;
        final int yN = y0 + targetRectangle.height - 1;
        final int x0 = targetRectangle.x;
        final int xN = targetRectangle.x + targetRectangle.width - 1;
        // System.out.println("x0 = " + x0 +", y0 = " + y0 + ", w = " + targetRectangle.width + ", h = " + targetRectangle.height);
        final int cohx0 = targetRectangle.x - (cohWinRg - 1) / 2;
        final int cohy0 = targetRectangle.y - (cohWinAz - 1) / 2;
        final int cohw = targetRectangle.width + cohWinRg - 1;
        final int cohh = targetRectangle.height + cohWinAz - 1;
        final Rectangle extRect = new Rectangle(cohx0, cohy0, cohw, cohh);
        final org.jlinda.core.Window tileWindow = new org.jlinda.core.Window(cohy0, cohy0 + cohh - 1, cohx0, cohx0 + cohw - 1);
        DemTile demTile = null;
        if (subtractTopographicPhase) {
            demTile = TopoPhase.getDEMTile(tileWindow, targetMap, dem, demNoDataValue, demSamplingLat, demSamplingLon, tileExtensionPercent);
            if (demTile.getData().length < 3 || demTile.getData()[0].length < 3) {
                throw new OperatorException("The resolution of the selected DEM is too low, " + "please select DEM with higher resolution.");
            }
        }
        for (String cohKey : targetMap.keySet()) {
            final ProductContainer product = targetMap.get(cohKey);
            final Tile mstTileReal = getSourceTile(product.sourceMaster.realBand, extRect, border);
            final Tile mstTileImag = getSourceTile(product.sourceMaster.imagBand, extRect, border);
            final ComplexDoubleMatrix dataMaster = TileUtilsDoris.pullComplexDoubleMatrix(mstTileReal, mstTileImag);
            final Tile slvTileReal = getSourceTile(product.sourceSlave.realBand, extRect, border);
            final Tile slvTileImag = getSourceTile(product.sourceSlave.imagBand, extRect, border);
            final ComplexDoubleMatrix dataSlave = TileUtilsDoris.pullComplexDoubleMatrix(slvTileReal, slvTileImag);
            if (subtractFlatEarthPhase) {
                final DoubleMatrix flatEarthPhase = computeFlatEarthPhase(cohx0, cohx0 + cohw - 1, cohw, cohy0, cohy0 + cohh - 1, cohh, 0, sourceImageWidth - 1, 0, sourceImageHeight - 1, product.sourceSlave.name);
                final ComplexDoubleMatrix complexReferencePhase = new ComplexDoubleMatrix(MatrixFunctions.cos(flatEarthPhase), MatrixFunctions.sin(flatEarthPhase));
                dataSlave.muli(complexReferencePhase);
                if (OUTPUT_PHASE) {
                    saveFlatEarthPhase(x0, xN, y0, yN, flatEarthPhase, product, targetTileMap);
                }
            }
            if (subtractTopographicPhase) {
                final TopoPhase topoPhase = TopoPhase.computeTopoPhase(product, tileWindow, demTile, false);
                final ComplexDoubleMatrix ComplexTopoPhase = new ComplexDoubleMatrix(MatrixFunctions.cos(new DoubleMatrix(topoPhase.demPhase)), MatrixFunctions.sin(new DoubleMatrix(topoPhase.demPhase)));
                dataSlave.muli(ComplexTopoPhase);
                if (OUTPUT_PHASE) {
                    saveTopoPhase(x0, xN, y0, yN, topoPhase.demPhase, product, targetTileMap);
                }
            }
            for (int i = 0; i < dataMaster.length; i++) {
                double tmp = norm(dataMaster.get(i));
                dataMaster.put(i, dataMaster.get(i).mul(dataSlave.get(i).conj()));
                dataSlave.put(i, new ComplexDouble(norm(dataSlave.get(i)), tmp));
            }
            DoubleMatrix cohMatrix = SarUtils.coherence2(dataMaster, dataSlave, cohWinAz, cohWinRg);
            saveCoherence(cohMatrix, product, targetTileMap, targetRectangle);
        }
    } catch (Throwable e) {
        OperatorUtils.catchOperatorException(getId(), e);
    } finally {
        pm.done();
    }
}
Also used : BorderExtender(javax.media.jai.BorderExtender) DemTile(org.jlinda.core.geom.DemTile) Tile(org.esa.snap.core.gpf.Tile) Point(org.jlinda.core.Point) GeoPoint(org.jlinda.core.GeoPoint) DoubleMatrix(org.jblas.DoubleMatrix) ComplexDoubleMatrix(org.jblas.ComplexDoubleMatrix) ComplexDouble(org.jblas.ComplexDouble) DemTile(org.jlinda.core.geom.DemTile) TopoPhase(org.jlinda.core.geom.TopoPhase) OperatorException(org.esa.snap.core.gpf.OperatorException) ComplexDoubleMatrix(org.jblas.ComplexDoubleMatrix)

Example 4 with DemTile

use of org.jlinda.core.geom.DemTile in project s1tbx by senbox-org.

the class InterferogramOp method computePartialTile.

private void computePartialTile(final int subSwathIndex, final int burstIndex, final int firstLineIdx, final Rectangle targetRectangle, final Map<Band, Tile> targetTileMap) {
    try {
        final BorderExtender border = BorderExtender.createInstance(BorderExtender.BORDER_ZERO);
        final int y0 = targetRectangle.y;
        final int yN = y0 + targetRectangle.height - 1;
        final int x0 = targetRectangle.x;
        final int xN = x0 + targetRectangle.width - 1;
        final Window tileWindow = new Window(y0 - firstLineIdx, yN - firstLineIdx, x0, xN);
        final SLCImage mstMeta = targetMap.values().iterator().next().sourceMaster.metaData.clone();
        updateMstMetaData(burstIndex, mstMeta);
        final Orbit mstOrbit = targetMap.values().iterator().next().sourceMaster.orbit;
        DemTile demTile = null;
        if (subtractTopographicPhase) {
            demTile = TopoPhase.getDEMTile(tileWindow, mstMeta, mstOrbit, dem, demNoDataValue, demSamplingLat, demSamplingLon, tileExtensionPercent);
            if (demTile == null) {
                throw new OperatorException("The selected DEM has no overlap with the image or is invalid.");
            }
            if (demTile.getData().length < 3 || demTile.getData()[0].length < 3) {
                throw new OperatorException("The resolution of the selected DEM is too low, " + "please select DEM with higher resolution.");
            }
        }
        final int cohx0 = targetRectangle.x - (cohWinRg - 1) / 2;
        final int cohy0 = targetRectangle.y - (cohWinAz - 1) / 2;
        final int cohw = targetRectangle.width + cohWinRg - 1;
        final int cohh = targetRectangle.height + cohWinAz - 1;
        final Rectangle rect = new Rectangle(cohx0, cohy0, cohw, cohh);
        final Window cohTileWindow = new Window(cohy0 - firstLineIdx, cohy0 + cohh - 1 - firstLineIdx, cohx0, cohx0 + cohw - 1);
        DemTile cohDemTile = null;
        if (subtractTopographicPhase) {
            cohDemTile = TopoPhase.getDEMTile(cohTileWindow, mstMeta, mstOrbit, dem, demNoDataValue, demSamplingLat, demSamplingLon, tileExtensionPercent);
        }
        final int minLine = 0;
        final int maxLine = subSwath[subSwathIndex - 1].linesPerBurst - 1;
        final int minPixel = 0;
        final int maxPixel = subSwath[subSwathIndex - 1].samplesPerBurst - 1;
        for (String ifgKey : targetMap.keySet()) {
            final ProductContainer product = targetMap.get(ifgKey);
            final SLCImage slvMeta = product.sourceSlave.metaData.clone();
            updateSlvMetaData(product, burstIndex, slvMeta);
            final Orbit slvOrbit = product.sourceSlave.orbit;
            // / check out results from master ///
            final Tile mstTileReal = getSourceTile(product.sourceMaster.realBand, targetRectangle, border);
            final Tile mstTileImag = getSourceTile(product.sourceMaster.imagBand, targetRectangle, border);
            final ComplexDoubleMatrix dataMaster = TileUtilsDoris.pullComplexDoubleMatrix(mstTileReal, mstTileImag);
            // / check out results from slave ///
            final Tile slvTileReal = getSourceTile(product.sourceSlave.realBand, targetRectangle, border);
            final Tile slvTileImag = getSourceTile(product.sourceSlave.imagBand, targetRectangle, border);
            final ComplexDoubleMatrix dataSlave = TileUtilsDoris.pullComplexDoubleMatrix(slvTileReal, slvTileImag);
            final String polynomialName = product.sourceSlave.name + '_' + (subSwathIndex - 1) + '_' + burstIndex;
            if (subtractFlatEarthPhase) {
                final DoubleMatrix flatEarthPhase = computeFlatEarthPhase(x0, xN, dataMaster.columns, y0 - firstLineIdx, yN - firstLineIdx, dataMaster.rows, minPixel, maxPixel, minLine, maxLine, polynomialName);
                final ComplexDoubleMatrix complexReferencePhase = new ComplexDoubleMatrix(MatrixFunctions.cos(flatEarthPhase), MatrixFunctions.sin(flatEarthPhase));
                dataSlave.muli(complexReferencePhase);
                if (OUTPUT_PHASE) {
                    saveFlatEarthPhase(x0, xN, y0, yN, flatEarthPhase, product, targetTileMap);
                }
            }
            if (subtractTopographicPhase) {
                TopoPhase topoPhase = TopoPhase.computeTopoPhase(mstMeta, mstOrbit, slvMeta, slvOrbit, tileWindow, demTile, outputElevation, false);
                final ComplexDoubleMatrix ComplexTopoPhase = new ComplexDoubleMatrix(MatrixFunctions.cos(new DoubleMatrix(topoPhase.demPhase)), MatrixFunctions.sin(new DoubleMatrix(topoPhase.demPhase)));
                dataSlave.muli(ComplexTopoPhase);
                if (OUTPUT_PHASE) {
                    saveTopoPhase(x0, xN, y0, yN, topoPhase.demPhase, product, targetTileMap);
                }
                if (outputElevation) {
                    saveElevation(x0, xN, y0, yN, topoPhase.elevation, product, targetTileMap);
                }
                if (outputLatLon) {
                    TopoPhase topoPhase1 = TopoPhase.computeTopoPhase(mstMeta, mstOrbit, slvMeta, slvOrbit, tileWindow, demTile, false, true);
                    saveLatLon(x0, xN, y0, yN, topoPhase1.latitude, topoPhase1.longitude, product, targetTileMap);
                }
            }
            dataMaster.muli(dataSlave.conji());
            saveInterferogram(dataMaster, product, targetTileMap, targetRectangle);
            // coherence calculation
            if (includeCoherence) {
                final Tile mstTileReal2 = getSourceTile(product.sourceMaster.realBand, rect, border);
                final Tile mstTileImag2 = getSourceTile(product.sourceMaster.imagBand, rect, border);
                final Tile slvTileReal2 = getSourceTile(product.sourceSlave.realBand, rect, border);
                final Tile slvTileImag2 = getSourceTile(product.sourceSlave.imagBand, rect, border);
                final ComplexDoubleMatrix dataMaster2 = TileUtilsDoris.pullComplexDoubleMatrix(mstTileReal2, mstTileImag2);
                final ComplexDoubleMatrix dataSlave2 = TileUtilsDoris.pullComplexDoubleMatrix(slvTileReal2, slvTileImag2);
                if (subtractFlatEarthPhase) {
                    final DoubleMatrix flatEarthPhase = computeFlatEarthPhase(cohx0, cohx0 + cohw - 1, cohw, cohy0 - firstLineIdx, cohy0 + cohh - 1 - firstLineIdx, cohh, minPixel, maxPixel, minLine, maxLine, polynomialName);
                    final ComplexDoubleMatrix complexReferencePhase = new ComplexDoubleMatrix(MatrixFunctions.cos(flatEarthPhase), MatrixFunctions.sin(flatEarthPhase));
                    dataSlave2.muli(complexReferencePhase);
                }
                if (subtractTopographicPhase) {
                    TopoPhase topoPhase = TopoPhase.computeTopoPhase(mstMeta, mstOrbit, slvMeta, slvOrbit, cohTileWindow, cohDemTile, false);
                    final ComplexDoubleMatrix ComplexTopoPhase = new ComplexDoubleMatrix(MatrixFunctions.cos(new DoubleMatrix(topoPhase.demPhase)), MatrixFunctions.sin(new DoubleMatrix(topoPhase.demPhase)));
                    dataSlave2.muli(ComplexTopoPhase);
                }
                for (int i = 0; i < dataMaster2.length; i++) {
                    double tmp = norm(dataMaster2.get(i));
                    dataMaster2.put(i, dataMaster2.get(i).mul(dataSlave2.get(i).conj()));
                    dataSlave2.put(i, new ComplexDouble(norm(dataSlave2.get(i)), tmp));
                }
                DoubleMatrix cohMatrix = SarUtils.coherence2(dataMaster2, dataSlave2, cohWinAz, cohWinRg);
                saveCoherence(cohMatrix, product, targetTileMap, targetRectangle);
            }
        }
    } catch (Throwable e) {
        OperatorUtils.catchOperatorException(getId(), e);
    }
}
Also used : Window(org.jlinda.core.Window) BorderExtender(javax.media.jai.BorderExtender) DemTile(org.jlinda.core.geom.DemTile) Tile(org.esa.snap.core.gpf.Tile) Point(org.jlinda.core.Point) DemTile(org.jlinda.core.geom.DemTile) TopoPhase(org.jlinda.core.geom.TopoPhase) OperatorException(org.esa.snap.core.gpf.OperatorException)

Example 5 with DemTile

use of org.jlinda.core.geom.DemTile in project s1tbx by senbox-org.

the class SimulateAmplitudeOp method computeTileStack.

/**
 * Called by the framework in order to compute a tile for the given target band.
 * <p>The default implementation throws a runtime exception with the message "not implemented".</p>
 *
 * @param targetTileMap   The target tiles associated with all target bands to be computed.
 * @param targetRectangle The rectangle of target tile.
 * @param pm              A progress monitor which should be used to determine computation cancelation requests.
 * @throws OperatorException
 *          If an error occurs during computation of the target raster.
 */
@Override
public void computeTileStack(Map<Band, Tile> targetTileMap, Rectangle targetRectangle, ProgressMonitor pm) throws OperatorException {
    try {
        int y0 = targetRectangle.y;
        int yN = y0 + targetRectangle.height - 1;
        int x0 = targetRectangle.x;
        int xN = targetRectangle.x + targetRectangle.width - 1;
        final Window tileWindow = new Window(y0, yN, x0, xN);
        Band simAmplitudeBand;
        Band targetBand_I;
        Band targetBand_Q;
        // need (extra space) for the coverage, taking into the consideration only height of the tile?
        for (String ifgKey : targetMap.keySet()) {
            ProductContainer product = targetMap.get(ifgKey);
            // / get dem of tile ///
            // compute tile geo-corners ~ work on ellipsoid
            GeoPoint[] geoCorners = GeoUtils.computeCorners(product.sourceMaster.metaData, product.sourceMaster.orbit, tileWindow);
            // get corners as DEM indices
            PixelPos[] pixelCorners = new PixelPos[2];
            pixelCorners[0] = dem.getIndex(new GeoPos(geoCorners[0].lat, geoCorners[0].lon));
            pixelCorners[1] = dem.getIndex(new GeoPos(geoCorners[1].lat, geoCorners[1].lon));
            // get max/min height of tile ~ uses 'fast' GCP based interpolation technique
            double[] tileHeights = computeMaxHeight(pixelCorners, targetRectangle);
            // compute extra lat/lon for dem tile
            GeoPoint geoExtent = GeoUtils.defineExtraPhiLam(tileHeights[0], tileHeights[1], tileWindow, product.sourceMaster.metaData, product.sourceMaster.orbit);
            // extend corners
            geoCorners = GeoUtils.extendCorners(geoExtent, geoCorners);
            // update corners
            pixelCorners[0] = dem.getIndex(new GeoPos(geoCorners[0].lat, geoCorners[0].lon));
            pixelCorners[1] = dem.getIndex(new GeoPos(geoCorners[1].lat, geoCorners[1].lon));
            pixelCorners[0] = new PixelPos(Math.ceil(pixelCorners[0].x), Math.floor(pixelCorners[0].y));
            pixelCorners[1] = new PixelPos(Math.floor(pixelCorners[1].x), Math.ceil(pixelCorners[1].y));
            GeoPos upperLeftGeo = dem.getGeoPos(pixelCorners[0]);
            int nLatPixels = (int) Math.abs(pixelCorners[1].y - pixelCorners[0].y);
            int nLonPixels = (int) Math.abs(pixelCorners[1].x - pixelCorners[0].x);
            int startX = (int) pixelCorners[0].x;
            int endX = startX + nLonPixels;
            int startY = (int) pixelCorners[0].y;
            int endY = startY + nLatPixels;
            double[][] elevation = new double[nLatPixels][nLonPixels];
            for (int y = startY, i = 0; y < endY; y++, i++) {
                for (int x = startX, j = 0; x < endX; x++, j++) {
                    try {
                        double elev = dem.getSample(x, y);
                        if (Double.isNaN(elev))
                            elev = demNoDataValue;
                        elevation[i][j] = elev;
                    } catch (Exception e) {
                        elevation[i][j] = demNoDataValue;
                    }
                }
            }
            DemTile demTile = new DemTile(upperLeftGeo.lat * Constants.DTOR, upperLeftGeo.lon * Constants.DTOR, nLatPixels, nLonPixels, demSampling, demSampling, (long) demNoDataValue);
            demTile.setData(elevation);
            final ThetaTile thetaTile = new ThetaTile(product.sourceMaster.metaData, product.sourceMaster.orbit, tileWindow, demTile);
            thetaTile.radarCode();
            thetaTile.gridData();
            final SimAmpTile simAmpTile = new SimAmpTile(product.sourceMaster.metaData, product.sourceMaster.orbit, tileWindow, thetaTile, demTile);
            simAmpTile.simulateAmplitude();
            // / check out results from source ///
            Tile tileReal = getSourceTile(product.sourceMaster.realBand, targetRectangle);
            Tile tileImag = getSourceTile(product.sourceMaster.imagBand, targetRectangle);
            ComplexDoubleMatrix complexIfg = TileUtilsDoris.pullComplexDoubleMatrix(tileReal, tileImag);
            // / commit to target ///
            targetBand_I = targetProduct.getBand(product.targetBandName_I);
            Tile tileOutReal = targetTileMap.get(targetBand_I);
            TileUtilsDoris.pushDoubleMatrix(complexIfg.real(), tileOutReal, targetRectangle);
            targetBand_Q = targetProduct.getBand(product.targetBandName_Q);
            Tile tileOutImag = targetTileMap.get(targetBand_Q);
            TileUtilsDoris.pushDoubleMatrix(complexIfg.imag(), tileOutImag, targetRectangle);
            // check in also simulated amplitude
            simAmplitudeBand = targetProduct.getBand(product.masterSubProduct.targetBandName_I);
            Tile tileOutSimAmplitude = targetTileMap.get(simAmplitudeBand);
            TileUtilsDoris.pushDoubleArray2D(simAmpTile.simAmpArray, tileOutSimAmplitude, targetRectangle);
        }
    } catch (Exception e) {
        throw new OperatorException(e);
    }
}
Also used : Window(org.jlinda.core.Window) DemTile(org.jlinda.core.geom.DemTile) Tile(org.esa.snap.core.gpf.Tile) SimAmpTile(org.jlinda.core.geom.SimAmpTile) ThetaTile(org.jlinda.core.geom.ThetaTile) OperatorException(org.esa.snap.core.gpf.OperatorException) IOException(java.io.IOException) DemTile(org.jlinda.core.geom.DemTile) SimAmpTile(org.jlinda.core.geom.SimAmpTile) ThetaTile(org.jlinda.core.geom.ThetaTile) OperatorException(org.esa.snap.core.gpf.OperatorException) ComplexDoubleMatrix(org.jblas.ComplexDoubleMatrix)

Aggregations

OperatorException (org.esa.snap.core.gpf.OperatorException)6 Tile (org.esa.snap.core.gpf.Tile)6 DemTile (org.jlinda.core.geom.DemTile)6 TopoPhase (org.jlinda.core.geom.TopoPhase)5 BorderExtender (javax.media.jai.BorderExtender)4 ComplexDoubleMatrix (org.jblas.ComplexDoubleMatrix)4 Point (org.jlinda.core.Point)4 Window (org.jlinda.core.Window)4 DoubleMatrix (org.jblas.DoubleMatrix)3 IOException (java.io.IOException)2 ComplexDouble (org.jblas.ComplexDouble)2 GeoPoint (org.jlinda.core.GeoPoint)2 Orbit (org.jlinda.core.Orbit)1 SLCImage (org.jlinda.core.SLCImage)1 SimAmpTile (org.jlinda.core.geom.SimAmpTile)1 ThetaTile (org.jlinda.core.geom.ThetaTile)1 ProductContainer (org.jlinda.core.utils.ProductContainer)1