Revision: 6625 http://sourceforge.net/p/jump-pilot/code/6625 Author: michaudm Date: 2020-11-19 08:02:57 +0000 (Thu, 19 Nov 2020) Log Message: ----------- Fix more deeply and document the half pixel shift problem
Modified Paths: -------------- core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoImage.java core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoReferencedRaster.java Modified: core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoImage.java =================================================================== --- core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoImage.java 2020-11-15 17:54:16 UTC (rev 6624) +++ core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoImage.java 2020-11-19 08:02:57 UTC (rev 6625) @@ -43,6 +43,7 @@ import java.awt.image.renderable.ParameterBlock; import java.lang.reflect.InvocationTargetException; +import javax.media.jai.InterpolationNearest; import javax.media.jai.JAI; import javax.media.jai.RenderedOp; import javax.media.jai.operator.AffineDescriptor; @@ -147,7 +148,7 @@ // System.out.println("GI: NO SCALE CACHE"); // First, scale the original image - double scaleX = scale * gtr.getDblModelUnitsPerRasterUnit_X(); + double scaleX = scale * Math.abs(gtr.getDblModelUnitsPerRasterUnit_X()); double scaleY = scale * Math.abs(gtr.getDblModelUnitsPerRasterUnit_Y()); // calculate predicted dimensions @@ -209,7 +210,11 @@ // Interpolation interp = Interpolation // .getInstance(Interpolation.INTERP_BICUBIC); // pb.add(interp); // add interpolation method - img = JAI.create("scale", pb, hints); + //img = JAI.create("scale", pb, hints); + // Prefer affine to scale as affine uses double parameters and preserve scale signus + // while scale uses float parameters and forbid negative scaleY + AffineTransform tr = new AffineTransform(scaleX_toUse, 0, 0, scaleY_toUse, 0, 0); + img = JAI.create("affine", scale_src_img, tr, new InterpolationNearest()); } else { pb.add(scaleX_toUse); pb.add(scaleY_toUse); @@ -270,9 +275,9 @@ raster_cropX = Math.min((float)raster_cropX, (float)img.getWidth()); raster_cropY = Math.min((float)raster_cropY, (float)img.getHeight()); raster_cropW = Math - .min((float)raster_cropW, (float)img.getWidth() - /*(int)*/ raster_cropX); - raster_cropH = Math.min((float)raster_cropH, (float)img.getHeight() - - /*(int)*/ raster_cropY); + .min((float)raster_cropW, (float)(img.getWidth() - raster_cropX)); + raster_cropH = Math.min((float)raster_cropH, (float)(img.getHeight() + - raster_cropY)); pb = new ParameterBlock(); pb.addSource(img); @@ -282,8 +287,11 @@ //System.out.println("croph " + (float)raster_cropH + " " + (img.getHeight() - /*(int)*/ raster_cropY)); pb.add((float)raster_cropX); pb.add((float)raster_cropY); - pb.add((float)raster_cropW); - pb.add((float)raster_cropH); + // conversions of cropX/Y and width/height from double to float may + // create an envelope slightly outside the real envelope of the image + // removing ulp to the width/height solve the problem + pb.add((float)raster_cropW-Math.ulp((float)raster_cropW)); + pb.add((float)raster_cropH-Math.ulp((float)raster_cropH)); img = JAI.create("crop", pb, null); // move the image to the model coordinates Modified: core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoReferencedRaster.java =================================================================== --- core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoReferencedRaster.java 2020-11-15 17:54:16 UTC (rev 6624) +++ core/trunk/src/com/vividsolutions/jump/workbench/imagery/geoimg/GeoReferencedRaster.java 2020-11-19 08:02:57 UTC (rev 6625) @@ -41,7 +41,6 @@ import java.io.IOException; import java.io.InputStream; import java.net.URI; -import java.net.URISyntaxException; import java.util.Arrays; import java.util.List; @@ -77,6 +76,10 @@ @Deprecated Coordinate coorModel_tiepointLT; + // AreaOrPoint=AREA (default) means that image coordinates + // refer to the upper left corner of the upper left angle. + // AreaOrPoint=POINT means that image coordinates refer to + // the center of the upper left pixel public enum AreaOrPoint {AREA,POINT} AreaOrPoint areaOrPoint = AreaOrPoint.AREA; @@ -91,6 +94,8 @@ private Coordinate modelULPixelCenter; private double dblModelUnitsPerRasterUnit_X; + // [michaudm 2020-11-16] signed y-scale + // to be able to handle north-south oriented or south-north oriented images private double dblModelUnitsPerRasterUnit_Y; /** @@ -120,7 +125,7 @@ re = new ReferencedImageException("probably no tiff image: " + e.getMessage()); } catch (IOException e) { - re = new ReferencedImageException("problem acessing tiff image: " + re = new ReferencedImageException("problem accessing tiff image: " + e.getMessage()); } finally { // clean up @@ -155,7 +160,7 @@ tags[4] = fieldModelGeoTransform.getAsDouble(3); // y-ordinate of the center of the upper left pixel tags[5] = fieldModelGeoTransform.getAsDouble(7); - Logger.debug("gtiff trans: " + Arrays.toString(tags)); + Logger.debug("gtiff transform: " + Arrays.toString(tags)); setEnvelope(tags); } // use the tiepoints as defined @@ -166,25 +171,20 @@ // Map the modeltiepoints from raster to model space // Read the tiepoints - double offset = (areaOrPoint == AreaOrPoint.AREA) ? 0.5 : 0.0; - // in areaOrPoint=AREA model, reference for the TiePoint is the UL - // corner of the pixel - // To set the transformation parameters we use the center of the - // UL coordinate + // imageCoord has not the same meaning in AREA and in POINT mode, + // but here, we don't mind, imageCoord may represent any point in + // the image, important thing is that imageCoord and modelCoord + // represent the same point in the image and in the model + // coordinate system. Coordinate imageCoord = new Coordinate( - fieldModelTiePoints.getAsDouble(0)-offset, - fieldModelTiePoints.getAsDouble(1)-offset); + fieldModelTiePoints.getAsDouble(0), + fieldModelTiePoints.getAsDouble(1)); Coordinate modelCoord = new Coordinate( fieldModelTiePoints.getAsDouble(3), fieldModelTiePoints.getAsDouble(4)); - //setCoorRasterTiff_tiepointLT(new Coordinate( - // fieldModelTiePoints.getAsDouble(0), - // fieldModelTiePoints.getAsDouble(1), 0)); - //setCoorModel_tiepointLT(new Coordinate( - // fieldModelTiePoints.getAsDouble(3), - // fieldModelTiePoints.getAsDouble(4), 0)); - //setEnvelope(); + Logger.debug("gtiff tiepoints found : " + Arrays.toString(fieldModelTiePoints.getAsDoubles())); + // Find the ModelPixelScale field XTIFFField fieldModelPixelScale = dir .getField(XTIFF.TIFFTAG_GEO_PIXEL_SCALE); @@ -195,8 +195,6 @@ + "\n" + MSG_GENERAL); } - Logger.debug("gtiff tiepoints found : " + Arrays.toString(fieldModelTiePoints.getAsDoubles())); - dblModelUnitsPerRasterUnit_X = fieldModelPixelScale.getAsDouble(0); dblModelUnitsPerRasterUnit_Y = fieldModelPixelScale.getAsDouble(1); //https://trac.osgeo.org/gdal/ticket/4977 @@ -203,14 +201,20 @@ if (!honourNegativeScaleY) dblModelUnitsPerRasterUnit_Y = - Math.abs(dblModelUnitsPerRasterUnit_Y); Logger.debug("gtiff scale : scalex=" + dblModelUnitsPerRasterUnit_X + ", scaley=" + dblModelUnitsPerRasterUnit_Y); - double tx = modelCoord.x-getDblModelUnitsPerRasterUnit_X()*(imageCoord.x); - double ty = modelCoord.y-getDblModelUnitsPerRasterUnit_Y()*(imageCoord.y); + // To compute the translation parameters of the transformation, we need + // to know how a point is converted from image to model coordinate, we + // don't mind where this point is. + double tx = modelCoord.x - dblModelUnitsPerRasterUnit_X * imageCoord.x; + double ty = modelCoord.y - dblModelUnitsPerRasterUnit_Y * imageCoord.y; - // Compute the model coordinate for the center of UL and LR pixels - rasterULPixelCenter = new Coordinate(0.0,0.0); + // Now, we want to know the model coordinate of the point precisely + // located at the center of the upper left pixel. + // Coordinates of this point is not the same in AREA and in POINT modes + rasterULPixelCenter = (areaOrPoint == AreaOrPoint.AREA) ? + new Coordinate(0.5,0.5) : new Coordinate(0.0,0.0); modelULPixelCenter = new Coordinate( - getDblModelUnitsPerRasterUnit_X()*rasterULPixelCenter.x + tx, - getDblModelUnitsPerRasterUnit_Y()*rasterULPixelCenter.y + ty); + getDblModelUnitsPerRasterUnit_X() * rasterULPixelCenter.x + tx, + getDblModelUnitsPerRasterUnit_Y() * rasterULPixelCenter.y + ty); setEnvelope(); @@ -294,7 +298,7 @@ tags[i] = Double.parseDouble(line); } - Logger.debug("wf: " + Arrays.toString(tags)); + Logger.debug("worldfile: " + Arrays.toString(tags)); setEnvelope(tags); @@ -355,24 +359,17 @@ // set up a default envelope double[] tags = new double[6]; - tags[0] = 1; // pixel size in x - // direction + tags[0] = 1; // pixel size in x direction tags[1] = 0; // rotation about y-axis tags[2] = 0; // rotation about x-axis - tags[3] = -1;// pixel size in the - // y-direction - tags[4] = 0; // x-coordinate of the - // center of the upper - // left pixel - tags[5] = 0; // y-coordinate of the - // center of the upper - // left pixel + tags[3] = -1;// pixel size in the y-direction + tags[4] = 0; // x-coordinate of the center of the upper left pixel + tags[5] = 0; // y-coordinate of the center of the upper left pixel setEnvelope(tags); } private void setEnvelope(double[] tags) { - //setCoorRasterTiff_tiepointLT(new Coordinate(-0.5, -0.5)); - //setCoorModel_tiepointLT(new Coordinate(0, 0)); + AffineTransform transform = new AffineTransform(tags); //We should honour negative scale y, but gdal created plenty of @@ -386,11 +383,14 @@ dblModelUnitsPerRasterUnit_X = transform.getScaleX(); dblModelUnitsPerRasterUnit_Y = transform.getScaleY(); - Point2D rasterLT = new Point2D.Double(src.getMinX(), src.getMinY()); + // To compute the envelope in AreaOrPoint.AREA mode, we need to + // know that upper left pixel center is at 0.5, 0.5, not 0, 0 + double offset = (areaOrPoint == AreaOrPoint.AREA) ? 0.5 : 0.0; + Point2D rasterLT = new Point2D.Double(src.getMinX()+offset, src.getMinY()+offset); Point2D modelLT = new Point2D.Double(); transform.transform(rasterLT, modelLT); - rasterULPixelCenter = new Coordinate(src.getMinX(), src.getMinY()); + rasterULPixelCenter = new Coordinate(rasterLT.getX(), rasterLT.getY()); modelULPixelCenter = new Coordinate(modelLT.getX(), modelLT.getY()); setEnvelope(); @@ -397,22 +397,23 @@ } void setEnvelope() { - double offset = (areaOrPoint == AreaOrPoint.AREA) ? 0.5 : 0.0; - // Get the image coordinate of the bottom left corner of the bottom left pixel - // from the image coordinate of the center of the bottom left pixel - Coordinate coorRaster_imageLB = new Coordinate( - rasterULPixelCenter.x-offset,src.getHeight()-offset); - // Get the image coordinate of the top right corner of the top right pixel - // from the image coordinate of the center of the top right pixel - Coordinate coorRaster_imageRT = new Coordinate( - src.getWidth()-offset, -offset); - // Transform the bottom left and the top right corner to model coordinate - // to get the envelope of the image - Coordinate coorModel_imageLB = rasterToModelSpace(coorRaster_imageLB); - Coordinate coorModel_imageRT = rasterToModelSpace(coorRaster_imageRT); - envModel_image = new Envelope(coorModel_imageLB, coorModel_imageRT); + // Image coordinate of the upper left corner of the envelope + double ulx = rasterULPixelCenter.x-0.5; + double uly = rasterULPixelCenter.y-0.5; + // Bottom left coordinate of the envelope + Coordinate imageEnvelopeBL = new Coordinate(ulx,uly+src.getHeight()); + + // Top right coordinate of the envelope + Coordinate imageEnvelopeTR = new Coordinate(ulx+src.getWidth(), uly); + + // Transform envelope corners to the model coordinate system + Coordinate modelEnvelopeBL = rasterToModelSpace(imageEnvelopeBL); + Coordinate modelEnvelopeTR = rasterToModelSpace(imageEnvelopeTR); + + envModel_image = new Envelope(modelEnvelopeBL, modelEnvelopeTR); + // backup original envelope envModel_image_backup = envModel_image; } @@ -426,13 +427,14 @@ */ private Coordinate rasterToModelSpace(Coordinate coorRaster) { Coordinate coorModel = new Coordinate(); - coorModel.x = modelULPixelCenter.x + (coorRaster.x - rasterULPixelCenter.x) * dblModelUnitsPerRasterUnit_X; + System.out.println("" + modelULPixelCenter.x + " + (" + coorRaster.x + "-" + rasterULPixelCenter.x + ")*" + dblModelUnitsPerRasterUnit_X + " = " + coorModel.x); coorModel.y = modelULPixelCenter.y + (coorRaster.y - rasterULPixelCenter.y) * dblModelUnitsPerRasterUnit_Y; + System.out.println("" + modelULPixelCenter.y + " + (" + coorRaster.y + "-" + rasterULPixelCenter.y + ")*" + dblModelUnitsPerRasterUnit_Y + " = " + coorModel.y); coorModel.z = 0; return coorModel; _______________________________________________ Jump-pilot-devel mailing list Jump-pilot-devel@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/jump-pilot-devel