@Test public void secondTest_fineMagFFT() throws Exception { processingPath = "/d2/test.processing/unit_tests/etna.volcano/process/crop/01486_21159.cpm/"; dataPath = "fineCoreg/magFFT/"; masterFileName = "master_patch.32x32.cr4"; slaveFileName = "slave_patch.32x32.cr4"; rows = 32; cols = 32; masterCplx = DataReader.readCplxFloatData(processingPath + dataPath + masterFileName, rows, cols, littleEndian); slaveCplx = DataReader.readCplxFloatData(processingPath + dataPath + slaveFileName, rows, cols, littleEndian); logger.trace("First entry of Master Complex Data Input {}: ", masterCplx.get(0, 0)); logger.trace("First entry of Slave Complex Data Input {}: ", slaveCplx.get(0, 0)); logger.trace("Master Input Dimensions [{},{}]", masterCplx.rows, masterCplx.columns); logger.trace("Slave Input Dimensions [{},{}]", slaveCplx.rows, slaveCplx.columns); clock.start(); Coregistration coreg = new Coregistration(); int MasksizeL = rows; int MasksizeP = cols; int ovsfactor; double coherence; for (int i = 0; i < 8; i++) { ovsfactor = (int) Math.pow(2, i); clock.start(); coherence = coreg.crosscorrelate(masterCplx, slaveCplx, ovsfactor, MasksizeL / 2, MasksizeP / 2, 0, 0); clock.stop(); logger.info("Estimated peak Coherence: {}, with OvsFactor: {}", coherence, ovsfactor); logger.info("Computation Time [ms]: {}", clock.getElapsedTime()); } }
@Override public double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL, int AccP, double offsetL, double offsetP) { logger.trace("crosscorrelate (PM 15-Apr-2012)"); final int L = Master.rows; final int P = Master.columns; final int twoL = 2 * L; final int twoP = 2 * P; final int halfL = L / 2; final int halfP = P / 2; if (Master.rows != Mask.rows || Master.columns != Mask.columns) { logger.error("Mask, Master not same size."); throw new IllegalArgumentException("Mask, Master not same size."); } if (!(MathUtils.isPower2(L) || MathUtils.isPower2(P))) { logger.error("Mask, Master size not power of 2."); throw new IllegalArgumentException("Mask, Master size not power of 2."); } if (!MathUtils.isPower2(ovsfactor)) { logger.error("coherencefft factor not power of 2"); throw new IllegalArgumentException("coherencefft factor not power of 2"); } logger.debug("Using de-meaned magnitude patches for incoherent cross-correlation"); DoubleMatrix magMaster = SarUtils.magnitude(Master); DoubleMatrix magMask = SarUtils.magnitude(Mask); magMaster.subi(magMaster.mean()); magMask.subi(magMask.mean()); ComplexDoubleMatrix Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); ComplexDoubleMatrix Mask2 = ComplexDoubleMatrix.zeros(twoL, twoP); Window windef = new Window(); Window win1 = new Window(0, L - 1, 0, P - 1); Window win2 = new Window(halfL, halfL + L - 1, halfP, halfP + P - 1); LinearAlgebraUtils.setdata(Master2, win1, new ComplexDoubleMatrix(magMaster), windef); LinearAlgebraUtils.setdata(Mask2, win2, new ComplexDoubleMatrix(magMask), windef); SpectralUtils.fft2D_inplace(Master2); SpectralUtils.fft2D_inplace(Mask2); Master2.conji(); Mask2.muli(Master2); SpectralUtils.invfft2D_inplace(Mask2); Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); int l, p; for (l = L; l < twoL; ++l) { for (p = P; p < twoP; ++p) { double realPart = magMaster.get(twoL - 1 - l, twoP - 1 - p); double imagPart = magMask.get(l - L, p - P); ComplexDouble value = new ComplexDouble(Math.pow(realPart, 2), Math.pow(imagPart, 2)); Master2.put(l, p, value); } } ComplexDoubleMatrix BLOCK = new ComplexDoubleMatrix(0, 0); if (BLOCK.rows != twoL || BLOCK.columns != twoP) { logger.debug("crosscorrelate:changing static block to size [" + twoL + ", " + twoP + "]"); BLOCK.resize(twoL, twoP); for (l = halfL; l < halfL + L; ++l) for (p = halfP; p < halfP + P; ++p) BLOCK.put(l, p, new ComplexDouble(1, 0)); SpectralUtils.fft2D_inplace(BLOCK); BLOCK.conji(); } SpectralUtils.fft2D_inplace(Master2); Master2.muli(BLOCK); SpectralUtils.invfft2D_inplace(Master2); DoubleMatrix Covar = new DoubleMatrix(L + 1, P + 1); double maxcorr = -999.0f; long maxcorrL = 0; long maxcorrP = 0; ComplexDouble maskValueTemp; ComplexDouble master2ValueTemp; for (l = 0; l <= L; ++l) { for (p = 0; p <= P; ++p) { maskValueTemp = Mask2.get(l, p); master2ValueTemp = Master2.get(l, p); Covar.put(l, p, maskValueTemp.real() / Math.sqrt(master2ValueTemp.real() * master2ValueTemp.imag())); if (Covar.get(l, p) > maxcorr) { maxcorr = Covar.get(l, p); maxcorrL = l; maxcorrP = p; } } } offsetL = -halfL + maxcorrL; offsetP = -halfP + maxcorrP; logger.debug("Pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); if (ovsfactor > 1) { if (maxcorrL < AccL) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = AccL; } if (maxcorrP < AccP) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = AccP; } if (maxcorrL > (L - AccL)) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = L - AccL; } if (maxcorrP > (P - AccP)) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = P - AccP; } Window win3 = new Window(maxcorrL - AccL, maxcorrL + AccL - 1, maxcorrP - AccP, maxcorrP + AccP - 1); final DoubleMatrix chip = new DoubleMatrix((int) win3.lines(), (int) win3.pixels()); LinearAlgebraUtils.setdata(chip, Covar, win3); int offL; int offP; DoubleMatrix chipOversampled = SarUtils.oversample(new ComplexDoubleMatrix(chip), ovsfactor, ovsfactor).getReal(); int corrIndex = chipOversampled.argmax(); offL = chipOversampled.indexColumns(corrIndex); offP = chipOversampled.indexRows(corrIndex); maxcorr = chipOversampled.get(corrIndex); offsetL = -halfL + maxcorrL - AccL + (double) offL / (double) ovsfactor; offsetP = -halfP + maxcorrP - AccP + (double) offP / (double) ovsfactor; System.out.println("Oversampling factor: " + ovsfactor); System.out.println("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); logger.debug("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); } return maxcorr; }
Coregistration implements ICoregistration { @Override public double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL, int AccP, double offsetL, double offsetP) { logger.trace("crosscorrelate (PM 15-Apr-2012)"); final int L = Master.rows; final int P = Master.columns; final int twoL = 2 * L; final int twoP = 2 * P; final int halfL = L / 2; final int halfP = P / 2; if (Master.rows != Mask.rows || Master.columns != Mask.columns) { logger.error("Mask, Master not same size."); throw new IllegalArgumentException("Mask, Master not same size."); } if (!(MathUtils.isPower2(L) || MathUtils.isPower2(P))) { logger.error("Mask, Master size not power of 2."); throw new IllegalArgumentException("Mask, Master size not power of 2."); } if (!MathUtils.isPower2(ovsfactor)) { logger.error("coherencefft factor not power of 2"); throw new IllegalArgumentException("coherencefft factor not power of 2"); } logger.debug("Using de-meaned magnitude patches for incoherent cross-correlation"); DoubleMatrix magMaster = SarUtils.magnitude(Master); DoubleMatrix magMask = SarUtils.magnitude(Mask); magMaster.subi(magMaster.mean()); magMask.subi(magMask.mean()); ComplexDoubleMatrix Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); ComplexDoubleMatrix Mask2 = ComplexDoubleMatrix.zeros(twoL, twoP); Window windef = new Window(); Window win1 = new Window(0, L - 1, 0, P - 1); Window win2 = new Window(halfL, halfL + L - 1, halfP, halfP + P - 1); LinearAlgebraUtils.setdata(Master2, win1, new ComplexDoubleMatrix(magMaster), windef); LinearAlgebraUtils.setdata(Mask2, win2, new ComplexDoubleMatrix(magMask), windef); SpectralUtils.fft2D_inplace(Master2); SpectralUtils.fft2D_inplace(Mask2); Master2.conji(); Mask2.muli(Master2); SpectralUtils.invfft2D_inplace(Mask2); Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); int l, p; for (l = L; l < twoL; ++l) { for (p = P; p < twoP; ++p) { double realPart = magMaster.get(twoL - 1 - l, twoP - 1 - p); double imagPart = magMask.get(l - L, p - P); ComplexDouble value = new ComplexDouble(Math.pow(realPart, 2), Math.pow(imagPart, 2)); Master2.put(l, p, value); } } ComplexDoubleMatrix BLOCK = new ComplexDoubleMatrix(0, 0); if (BLOCK.rows != twoL || BLOCK.columns != twoP) { logger.debug("crosscorrelate:changing static block to size [" + twoL + ", " + twoP + "]"); BLOCK.resize(twoL, twoP); for (l = halfL; l < halfL + L; ++l) for (p = halfP; p < halfP + P; ++p) BLOCK.put(l, p, new ComplexDouble(1, 0)); SpectralUtils.fft2D_inplace(BLOCK); BLOCK.conji(); } SpectralUtils.fft2D_inplace(Master2); Master2.muli(BLOCK); SpectralUtils.invfft2D_inplace(Master2); DoubleMatrix Covar = new DoubleMatrix(L + 1, P + 1); double maxcorr = -999.0f; long maxcorrL = 0; long maxcorrP = 0; ComplexDouble maskValueTemp; ComplexDouble master2ValueTemp; for (l = 0; l <= L; ++l) { for (p = 0; p <= P; ++p) { maskValueTemp = Mask2.get(l, p); master2ValueTemp = Master2.get(l, p); Covar.put(l, p, maskValueTemp.real() / Math.sqrt(master2ValueTemp.real() * master2ValueTemp.imag())); if (Covar.get(l, p) > maxcorr) { maxcorr = Covar.get(l, p); maxcorrL = l; maxcorrP = p; } } } offsetL = -halfL + maxcorrL; offsetP = -halfP + maxcorrP; logger.debug("Pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); if (ovsfactor > 1) { if (maxcorrL < AccL) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = AccL; } if (maxcorrP < AccP) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = AccP; } if (maxcorrL > (L - AccL)) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = L - AccL; } if (maxcorrP > (P - AccP)) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = P - AccP; } Window win3 = new Window(maxcorrL - AccL, maxcorrL + AccL - 1, maxcorrP - AccP, maxcorrP + AccP - 1); final DoubleMatrix chip = new DoubleMatrix((int) win3.lines(), (int) win3.pixels()); LinearAlgebraUtils.setdata(chip, Covar, win3); int offL; int offP; DoubleMatrix chipOversampled = SarUtils.oversample(new ComplexDoubleMatrix(chip), ovsfactor, ovsfactor).getReal(); int corrIndex = chipOversampled.argmax(); offL = chipOversampled.indexColumns(corrIndex); offP = chipOversampled.indexRows(corrIndex); maxcorr = chipOversampled.get(corrIndex); offsetL = -halfL + maxcorrL - AccL + (double) offL / (double) ovsfactor; offsetP = -halfP + maxcorrP - AccP + (double) offP / (double) ovsfactor; System.out.println("Oversampling factor: " + ovsfactor); System.out.println("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); logger.debug("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); } return maxcorr; } }
Coregistration implements ICoregistration { @Override public double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL, int AccP, double offsetL, double offsetP) { logger.trace("crosscorrelate (PM 15-Apr-2012)"); final int L = Master.rows; final int P = Master.columns; final int twoL = 2 * L; final int twoP = 2 * P; final int halfL = L / 2; final int halfP = P / 2; if (Master.rows != Mask.rows || Master.columns != Mask.columns) { logger.error("Mask, Master not same size."); throw new IllegalArgumentException("Mask, Master not same size."); } if (!(MathUtils.isPower2(L) || MathUtils.isPower2(P))) { logger.error("Mask, Master size not power of 2."); throw new IllegalArgumentException("Mask, Master size not power of 2."); } if (!MathUtils.isPower2(ovsfactor)) { logger.error("coherencefft factor not power of 2"); throw new IllegalArgumentException("coherencefft factor not power of 2"); } logger.debug("Using de-meaned magnitude patches for incoherent cross-correlation"); DoubleMatrix magMaster = SarUtils.magnitude(Master); DoubleMatrix magMask = SarUtils.magnitude(Mask); magMaster.subi(magMaster.mean()); magMask.subi(magMask.mean()); ComplexDoubleMatrix Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); ComplexDoubleMatrix Mask2 = ComplexDoubleMatrix.zeros(twoL, twoP); Window windef = new Window(); Window win1 = new Window(0, L - 1, 0, P - 1); Window win2 = new Window(halfL, halfL + L - 1, halfP, halfP + P - 1); LinearAlgebraUtils.setdata(Master2, win1, new ComplexDoubleMatrix(magMaster), windef); LinearAlgebraUtils.setdata(Mask2, win2, new ComplexDoubleMatrix(magMask), windef); SpectralUtils.fft2D_inplace(Master2); SpectralUtils.fft2D_inplace(Mask2); Master2.conji(); Mask2.muli(Master2); SpectralUtils.invfft2D_inplace(Mask2); Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); int l, p; for (l = L; l < twoL; ++l) { for (p = P; p < twoP; ++p) { double realPart = magMaster.get(twoL - 1 - l, twoP - 1 - p); double imagPart = magMask.get(l - L, p - P); ComplexDouble value = new ComplexDouble(Math.pow(realPart, 2), Math.pow(imagPart, 2)); Master2.put(l, p, value); } } ComplexDoubleMatrix BLOCK = new ComplexDoubleMatrix(0, 0); if (BLOCK.rows != twoL || BLOCK.columns != twoP) { logger.debug("crosscorrelate:changing static block to size [" + twoL + ", " + twoP + "]"); BLOCK.resize(twoL, twoP); for (l = halfL; l < halfL + L; ++l) for (p = halfP; p < halfP + P; ++p) BLOCK.put(l, p, new ComplexDouble(1, 0)); SpectralUtils.fft2D_inplace(BLOCK); BLOCK.conji(); } SpectralUtils.fft2D_inplace(Master2); Master2.muli(BLOCK); SpectralUtils.invfft2D_inplace(Master2); DoubleMatrix Covar = new DoubleMatrix(L + 1, P + 1); double maxcorr = -999.0f; long maxcorrL = 0; long maxcorrP = 0; ComplexDouble maskValueTemp; ComplexDouble master2ValueTemp; for (l = 0; l <= L; ++l) { for (p = 0; p <= P; ++p) { maskValueTemp = Mask2.get(l, p); master2ValueTemp = Master2.get(l, p); Covar.put(l, p, maskValueTemp.real() / Math.sqrt(master2ValueTemp.real() * master2ValueTemp.imag())); if (Covar.get(l, p) > maxcorr) { maxcorr = Covar.get(l, p); maxcorrL = l; maxcorrP = p; } } } offsetL = -halfL + maxcorrL; offsetP = -halfP + maxcorrP; logger.debug("Pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); if (ovsfactor > 1) { if (maxcorrL < AccL) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = AccL; } if (maxcorrP < AccP) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = AccP; } if (maxcorrL > (L - AccL)) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = L - AccL; } if (maxcorrP > (P - AccP)) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = P - AccP; } Window win3 = new Window(maxcorrL - AccL, maxcorrL + AccL - 1, maxcorrP - AccP, maxcorrP + AccP - 1); final DoubleMatrix chip = new DoubleMatrix((int) win3.lines(), (int) win3.pixels()); LinearAlgebraUtils.setdata(chip, Covar, win3); int offL; int offP; DoubleMatrix chipOversampled = SarUtils.oversample(new ComplexDoubleMatrix(chip), ovsfactor, ovsfactor).getReal(); int corrIndex = chipOversampled.argmax(); offL = chipOversampled.indexColumns(corrIndex); offP = chipOversampled.indexRows(corrIndex); maxcorr = chipOversampled.get(corrIndex); offsetL = -halfL + maxcorrL - AccL + (double) offL / (double) ovsfactor; offsetP = -halfP + maxcorrP - AccP + (double) offP / (double) ovsfactor; System.out.println("Oversampling factor: " + ovsfactor); System.out.println("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); logger.debug("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); } return maxcorr; }  }
Coregistration implements ICoregistration { @Override public double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL, int AccP, double offsetL, double offsetP) { logger.trace("crosscorrelate (PM 15-Apr-2012)"); final int L = Master.rows; final int P = Master.columns; final int twoL = 2 * L; final int twoP = 2 * P; final int halfL = L / 2; final int halfP = P / 2; if (Master.rows != Mask.rows || Master.columns != Mask.columns) { logger.error("Mask, Master not same size."); throw new IllegalArgumentException("Mask, Master not same size."); } if (!(MathUtils.isPower2(L) || MathUtils.isPower2(P))) { logger.error("Mask, Master size not power of 2."); throw new IllegalArgumentException("Mask, Master size not power of 2."); } if (!MathUtils.isPower2(ovsfactor)) { logger.error("coherencefft factor not power of 2"); throw new IllegalArgumentException("coherencefft factor not power of 2"); } logger.debug("Using de-meaned magnitude patches for incoherent cross-correlation"); DoubleMatrix magMaster = SarUtils.magnitude(Master); DoubleMatrix magMask = SarUtils.magnitude(Mask); magMaster.subi(magMaster.mean()); magMask.subi(magMask.mean()); ComplexDoubleMatrix Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); ComplexDoubleMatrix Mask2 = ComplexDoubleMatrix.zeros(twoL, twoP); Window windef = new Window(); Window win1 = new Window(0, L - 1, 0, P - 1); Window win2 = new Window(halfL, halfL + L - 1, halfP, halfP + P - 1); LinearAlgebraUtils.setdata(Master2, win1, new ComplexDoubleMatrix(magMaster), windef); LinearAlgebraUtils.setdata(Mask2, win2, new ComplexDoubleMatrix(magMask), windef); SpectralUtils.fft2D_inplace(Master2); SpectralUtils.fft2D_inplace(Mask2); Master2.conji(); Mask2.muli(Master2); SpectralUtils.invfft2D_inplace(Mask2); Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); int l, p; for (l = L; l < twoL; ++l) { for (p = P; p < twoP; ++p) { double realPart = magMaster.get(twoL - 1 - l, twoP - 1 - p); double imagPart = magMask.get(l - L, p - P); ComplexDouble value = new ComplexDouble(Math.pow(realPart, 2), Math.pow(imagPart, 2)); Master2.put(l, p, value); } } ComplexDoubleMatrix BLOCK = new ComplexDoubleMatrix(0, 0); if (BLOCK.rows != twoL || BLOCK.columns != twoP) { logger.debug("crosscorrelate:changing static block to size [" + twoL + ", " + twoP + "]"); BLOCK.resize(twoL, twoP); for (l = halfL; l < halfL + L; ++l) for (p = halfP; p < halfP + P; ++p) BLOCK.put(l, p, new ComplexDouble(1, 0)); SpectralUtils.fft2D_inplace(BLOCK); BLOCK.conji(); } SpectralUtils.fft2D_inplace(Master2); Master2.muli(BLOCK); SpectralUtils.invfft2D_inplace(Master2); DoubleMatrix Covar = new DoubleMatrix(L + 1, P + 1); double maxcorr = -999.0f; long maxcorrL = 0; long maxcorrP = 0; ComplexDouble maskValueTemp; ComplexDouble master2ValueTemp; for (l = 0; l <= L; ++l) { for (p = 0; p <= P; ++p) { maskValueTemp = Mask2.get(l, p); master2ValueTemp = Master2.get(l, p); Covar.put(l, p, maskValueTemp.real() / Math.sqrt(master2ValueTemp.real() * master2ValueTemp.imag())); if (Covar.get(l, p) > maxcorr) { maxcorr = Covar.get(l, p); maxcorrL = l; maxcorrP = p; } } } offsetL = -halfL + maxcorrL; offsetP = -halfP + maxcorrP; logger.debug("Pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); if (ovsfactor > 1) { if (maxcorrL < AccL) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = AccL; } if (maxcorrP < AccP) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = AccP; } if (maxcorrL > (L - AccL)) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = L - AccL; } if (maxcorrP > (P - AccP)) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = P - AccP; } Window win3 = new Window(maxcorrL - AccL, maxcorrL + AccL - 1, maxcorrP - AccP, maxcorrP + AccP - 1); final DoubleMatrix chip = new DoubleMatrix((int) win3.lines(), (int) win3.pixels()); LinearAlgebraUtils.setdata(chip, Covar, win3); int offL; int offP; DoubleMatrix chipOversampled = SarUtils.oversample(new ComplexDoubleMatrix(chip), ovsfactor, ovsfactor).getReal(); int corrIndex = chipOversampled.argmax(); offL = chipOversampled.indexColumns(corrIndex); offP = chipOversampled.indexRows(corrIndex); maxcorr = chipOversampled.get(corrIndex); offsetL = -halfL + maxcorrL - AccL + (double) offL / (double) ovsfactor; offsetP = -halfP + maxcorrP - AccP + (double) offP / (double) ovsfactor; System.out.println("Oversampling factor: " + ovsfactor); System.out.println("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); logger.debug("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); } return maxcorr; }  @Override void coarseporbit(Ellipsoid ell, SLCImage master, SLCImage slave, Orbit masterorbit, Orbit slaveorbit, Baseline baseline); @Override void coarsecorrel(Input.CoarseCorr input, SLCImage minfo, SLCImage sinfo); @Override void coarsecorrelfft(Input.CoarseCorr coarsecorrinput, SLCImage minfo, SLCImage sinfo); @Override void mtiming_correl(Input.MasterTiming input, SLCImage minfo, todo_classes.productinfo sinfo); @Override void mtiming_correlfft(Input.MasterTiming input, SLCImage minfo, todo_classes.productinfo sinfo); @Override int[][] distributepoints(float numberofpoints, Window win); @Override void getoffset(DoubleMatrix Result, long offsetLines, long offsetPixels); @Override void getmodeoffset(float[][] Result, long offsetLines, long offsetPixels); @Override void finecoreg(Input.FineCorr fineinput, SLCImage minfo, SLCImage sinfo); @Override double coherencefft(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL,
                               int AccP, double offsetL, double offsetP); @Override double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask,
                                 int ovsfactor,
                                 int AccL, int AccP,
                                 double offsetL, double offsetP); @Override double coherencespace(final int AccL, final int AccP, final int osfactor,
                                 ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask,
                                 double offsetL, double offsetP); @Override void coregpm(SLCImage master, SLCImage slave, String i_resfile, Input.CoregPM coregpminput,
                        int demassist); @Override DoubleMatrix getofffile(String file, float threshold); @Override void resample(Input.Resample resampleinput,
                         SLCImage master, SLCImage slave,
                         double[] cpmL, double[] cpmP,
                         int demassist); @Override void ms_timing_error(SLCImage master, String i_resfile, Input.RelTiming timinginput,
                                int coarse_orbit_offsetL, int coarse_orbit_offsetP); DoubleMatrix correlate(DoubleMatrix A, DoubleMatrix Mask); void shiftazispectrum(final ComplexDoubleMatrix data, // slcData in space domain
                                 final double prf,
                                 final double rsr2x,
                                 final double[] doppler,   // polynomial fDC
                                 final double shift); @Deprecated static Window getOverlap(final Window master, final Window slave); static Window getOverlap(final SLCImage master, final SLCImage slave, final double[] cpmL, final double[] cpmP); static Window getOverlap(final SLCImage master, final SLCImage slave,
                                    final double Npointsd2, final double timing_L, final double timing_P); }
Coregistration implements ICoregistration { @Override public double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL, int AccP, double offsetL, double offsetP) { logger.trace("crosscorrelate (PM 15-Apr-2012)"); final int L = Master.rows; final int P = Master.columns; final int twoL = 2 * L; final int twoP = 2 * P; final int halfL = L / 2; final int halfP = P / 2; if (Master.rows != Mask.rows || Master.columns != Mask.columns) { logger.error("Mask, Master not same size."); throw new IllegalArgumentException("Mask, Master not same size."); } if (!(MathUtils.isPower2(L) || MathUtils.isPower2(P))) { logger.error("Mask, Master size not power of 2."); throw new IllegalArgumentException("Mask, Master size not power of 2."); } if (!MathUtils.isPower2(ovsfactor)) { logger.error("coherencefft factor not power of 2"); throw new IllegalArgumentException("coherencefft factor not power of 2"); } logger.debug("Using de-meaned magnitude patches for incoherent cross-correlation"); DoubleMatrix magMaster = SarUtils.magnitude(Master); DoubleMatrix magMask = SarUtils.magnitude(Mask); magMaster.subi(magMaster.mean()); magMask.subi(magMask.mean()); ComplexDoubleMatrix Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); ComplexDoubleMatrix Mask2 = ComplexDoubleMatrix.zeros(twoL, twoP); Window windef = new Window(); Window win1 = new Window(0, L - 1, 0, P - 1); Window win2 = new Window(halfL, halfL + L - 1, halfP, halfP + P - 1); LinearAlgebraUtils.setdata(Master2, win1, new ComplexDoubleMatrix(magMaster), windef); LinearAlgebraUtils.setdata(Mask2, win2, new ComplexDoubleMatrix(magMask), windef); SpectralUtils.fft2D_inplace(Master2); SpectralUtils.fft2D_inplace(Mask2); Master2.conji(); Mask2.muli(Master2); SpectralUtils.invfft2D_inplace(Mask2); Master2 = ComplexDoubleMatrix.zeros(twoL, twoP); int l, p; for (l = L; l < twoL; ++l) { for (p = P; p < twoP; ++p) { double realPart = magMaster.get(twoL - 1 - l, twoP - 1 - p); double imagPart = magMask.get(l - L, p - P); ComplexDouble value = new ComplexDouble(Math.pow(realPart, 2), Math.pow(imagPart, 2)); Master2.put(l, p, value); } } ComplexDoubleMatrix BLOCK = new ComplexDoubleMatrix(0, 0); if (BLOCK.rows != twoL || BLOCK.columns != twoP) { logger.debug("crosscorrelate:changing static block to size [" + twoL + ", " + twoP + "]"); BLOCK.resize(twoL, twoP); for (l = halfL; l < halfL + L; ++l) for (p = halfP; p < halfP + P; ++p) BLOCK.put(l, p, new ComplexDouble(1, 0)); SpectralUtils.fft2D_inplace(BLOCK); BLOCK.conji(); } SpectralUtils.fft2D_inplace(Master2); Master2.muli(BLOCK); SpectralUtils.invfft2D_inplace(Master2); DoubleMatrix Covar = new DoubleMatrix(L + 1, P + 1); double maxcorr = -999.0f; long maxcorrL = 0; long maxcorrP = 0; ComplexDouble maskValueTemp; ComplexDouble master2ValueTemp; for (l = 0; l <= L; ++l) { for (p = 0; p <= P; ++p) { maskValueTemp = Mask2.get(l, p); master2ValueTemp = Master2.get(l, p); Covar.put(l, p, maskValueTemp.real() / Math.sqrt(master2ValueTemp.real() * master2ValueTemp.imag())); if (Covar.get(l, p) > maxcorr) { maxcorr = Covar.get(l, p); maxcorrL = l; maxcorrP = p; } } } offsetL = -halfL + maxcorrL; offsetP = -halfP + maxcorrP; logger.debug("Pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); if (ovsfactor > 1) { if (maxcorrL < AccL) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = AccL; } if (maxcorrP < AccP) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = AccP; } if (maxcorrL > (L - AccL)) { logger.debug("Careful, decrease AccL or increase winsizeL"); maxcorrL = L - AccL; } if (maxcorrP > (P - AccP)) { logger.debug("Careful, decrease AccP or increase winsizeP"); maxcorrP = P - AccP; } Window win3 = new Window(maxcorrL - AccL, maxcorrL + AccL - 1, maxcorrP - AccP, maxcorrP + AccP - 1); final DoubleMatrix chip = new DoubleMatrix((int) win3.lines(), (int) win3.pixels()); LinearAlgebraUtils.setdata(chip, Covar, win3); int offL; int offP; DoubleMatrix chipOversampled = SarUtils.oversample(new ComplexDoubleMatrix(chip), ovsfactor, ovsfactor).getReal(); int corrIndex = chipOversampled.argmax(); offL = chipOversampled.indexColumns(corrIndex); offP = chipOversampled.indexRows(corrIndex); maxcorr = chipOversampled.get(corrIndex); offsetL = -halfL + maxcorrL - AccL + (double) offL / (double) ovsfactor; offsetP = -halfP + maxcorrP - AccP + (double) offP / (double) ovsfactor; System.out.println("Oversampling factor: " + ovsfactor); System.out.println("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); logger.debug("Sub-pixel level offset: " + offsetL + ", " + offsetP + " (corr=" + maxcorr + ")"); } return maxcorr; }  @Override void coarseporbit(Ellipsoid ell, SLCImage master, SLCImage slave, Orbit masterorbit, Orbit slaveorbit, Baseline baseline); @Override void coarsecorrel(Input.CoarseCorr input, SLCImage minfo, SLCImage sinfo); @Override void coarsecorrelfft(Input.CoarseCorr coarsecorrinput, SLCImage minfo, SLCImage sinfo); @Override void mtiming_correl(Input.MasterTiming input, SLCImage minfo, todo_classes.productinfo sinfo); @Override void mtiming_correlfft(Input.MasterTiming input, SLCImage minfo, todo_classes.productinfo sinfo); @Override int[][] distributepoints(float numberofpoints, Window win); @Override void getoffset(DoubleMatrix Result, long offsetLines, long offsetPixels); @Override void getmodeoffset(float[][] Result, long offsetLines, long offsetPixels); @Override void finecoreg(Input.FineCorr fineinput, SLCImage minfo, SLCImage sinfo); @Override double coherencefft(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask, int ovsfactor, int AccL,
                               int AccP, double offsetL, double offsetP); @Override double crosscorrelate(ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask,
                                 int ovsfactor,
                                 int AccL, int AccP,
                                 double offsetL, double offsetP); @Override double coherencespace(final int AccL, final int AccP, final int osfactor,
                                 ComplexDoubleMatrix Master, ComplexDoubleMatrix Mask,
                                 double offsetL, double offsetP); @Override void coregpm(SLCImage master, SLCImage slave, String i_resfile, Input.CoregPM coregpminput,
                        int demassist); @Override DoubleMatrix getofffile(String file, float threshold); @Override void resample(Input.Resample resampleinput,
                         SLCImage master, SLCImage slave,
                         double[] cpmL, double[] cpmP,
                         int demassist); @Override void ms_timing_error(SLCImage master, String i_resfile, Input.RelTiming timinginput,
                                int coarse_orbit_offsetL, int coarse_orbit_offsetP); DoubleMatrix correlate(DoubleMatrix A, DoubleMatrix Mask); void shiftazispectrum(final ComplexDoubleMatrix data, // slcData in space domain
                                 final double prf,
                                 final double rsr2x,
                                 final double[] doppler,   // polynomial fDC
                                 final double shift); @Deprecated static Window getOverlap(final Window master, final Window slave); static Window getOverlap(final SLCImage master, final SLCImage slave, final double[] cpmL, final double[] cpmP); static Window getOverlap(final SLCImage master, final SLCImage slave,
                                    final double Npointsd2, final double timing_L, final double timing_P);  }