Mercurial > dive4elements > gnv-client
changeset 547:23d5cc37dd5b
Fixed access to raster data.
geo-backend/trunk@518 c6561f87-3c4e-4783-a992-168aeb5c3f6f
| author | Sascha L. Teichmann <sascha.teichmann@intevation.de> |
|---|---|
| date | Sat, 09 Jan 2010 12:41:26 +0000 |
| parents | 210716612c30 |
| children | ccd976fc0f7b |
| files | geo-backend/ChangeLog geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/ArcSDEStatement.java geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java |
| diffstat | 3 files changed, 231 insertions(+), 123 deletions(-) [+] |
line wrap: on
line diff
--- a/geo-backend/ChangeLog Fri Jan 08 14:40:08 2010 +0000 +++ b/geo-backend/ChangeLog Sat Jan 09 12:41:26 2010 +0000 @@ -1,3 +1,18 @@ +2009-01-09 Sascha L. Teichmann <sascha.teichmann@intevation.de> + + * src/main/java/de/intevation/gnv/geobackend/sde/datasources/ArcSDEStatement.java: + Fixed transformations to figure out the right tile in the SDE. + Taking the tile origin into account now and does not do wrong rounding + of column and row indices. Simplified the code a lot. + + * src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java: + Fixed access to right right pixel now. Removed the envelope code because + its not needed any longer. You can ask a RasterObject directly now if + a given coordinate is on the tile now. The tile space is now correctly + spanned between (0, 0) and (tile width - 1, tile height - 1) now which + allows other access patterns like bilinear interpolation. TODO: implement + bilinear interpolation. + 2010-01-08 Tim Englich <tim.englich@intevation.de> * src/main/java/de/intevation/gnv/geobackend/sde/datasources/SDEQuery.java,
--- a/geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/ArcSDEStatement.java Fri Jan 08 14:40:08 2010 +0000 +++ b/geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/ArcSDEStatement.java Sat Jan 09 12:41:26 2010 +0000 @@ -37,9 +37,9 @@ import com.vividsolutions.jts.io.WKTReader; - /** - * @author Tim Englich <tim.englich@intevation.de> - * +/** + * @author Tim Englich (tim.englich@intevation.de) + * @author Sascha L. Teichmann (sascha.teichmann@intevation.de) */ public class ArcSDEStatement implements Statement { @@ -484,7 +484,13 @@ * @throws SeException */ private ResultSet handleResultSet(SeQuery pSeQuery, boolean isRaster, Geometry geometry) throws SeException { - log.debug("ArcSDEStatement,handleResultSet()"); + + boolean debug = log.isDebugEnabled(); + + if (debug) { + log.debug("ArcSDEStatement.handleResultSet()"); + } + SDEResultSet lSet = new SDEResultSet(); SeRow row = null;; int lCount; @@ -495,7 +501,8 @@ // analyze cols of result set SeColumnDefinition[] lCols = row.getColumns(); for (SeColumnDefinition lCol : lCols) { - lSet.addCol(new ColDefinition(lCol.getName(), lCol.getType()));// notice: esri-types have been copied into colDefinition class! + lSet.addCol(new ColDefinition(lCol.getName(), lCol.getType())); + // notice: esri-types have been copied into colDefinition class! } } short lNumCols = row.getNumColumns(); @@ -507,44 +514,30 @@ } }else{ try { - - int rasterSize = 128; pSeQuery.execute(); row = pSeQuery.fetch(); - SeRasterAttr attr = row.getRaster(0); - SeRaster raster = attr.getRasterInfo(); - SeRasterBand[] bands = raster.getBands(); - SeRasterBand rasterBand = bands[0]; - - - int height = rasterBand.getBandHeight(); - int width = rasterBand.getBandWidth(); + SeRasterAttr attr = row.getRaster(0); + SeRaster raster = attr.getRasterInfo(); + SeRasterBand [] bands = raster.getBands(); + SeRasterBand rasterBand = bands[0]; SeExtent extent = rasterBand.getExtent(); - log.debug("Extent: "+ extent.getMinX() +" "+extent.getMinY()+" "+extent.getMaxX()+ " "+extent.getMaxY() ); - log.debug("Querygeometry: "+geometry.toText()); + + /* + if (debug) { + log.debug("Extent: " + + extent.getMinX() + " " + extent.getMinY() + " "+ + extent.getMaxX() + " " + extent.getMaxY()); + log.debug("Query geometry: "+geometry.toText()); + } + */ double x = ((Point)geometry).getX(); double y = ((Point)geometry).getY(); - boolean isPointInRaster = false; - if ((x >= extent.getMinX() && x <= extent.getMaxX()) && (y >= extent.getMinY() && y <= extent.getMaxY())){ - isPointInRaster = true; - } - - double dxNature = x - extent.getMinX(); - double dyNature = extent.getMaxY() - y;//y-extent.getMinY(); - double widthNature = extent.getMaxX()-extent.getMinX(); - double heightNature = extent.getMaxY()-extent.getMinY(); - - int pixelX = (int)Math.round(dxNature * width / widthNature); - int pixelY = (int)Math.round(dyNature * height/ heightNature); - - log.debug("Rasterdimesion (Pixel) " + width + " / " + height); - log.debug ("Required Pixel: " + pixelX + " / "+ pixelY); - int bandNumber = rasterBand.getBandNumber(); - log.debug("BAND-ID "+rasterBand.getId().longValue()); - int maxLevel = 0; - log.debug("Maxlevel: "+maxLevel); + + boolean isPointInRaster = + x >= extent.getMinX() && x <= extent.getMaxX() + && y >= extent.getMinY() && y <= extent.getMaxY(); if (isPointInRaster){ @@ -552,41 +545,111 @@ pSeQuery.execute(); row = pSeQuery.fetch(); } + + double midX = 0.5d*(extent.getMinX() + extent.getMaxX()); + double midY = 0.5d*(extent.getMinY() + extent.getMaxY()); + + SDEPoint origin = rasterBand.getTileOrigin(); + + double maxX = origin.getX() < midX + ? extent.getMaxX() + : extent.getMinX(); + + double maxY = origin.getY() < midY + ? extent.getMaxY() + : extent.getMinY(); + + /* + 0 = origin.getX()*mx + bx + rasterBand.getBandWidth() = maxX + bx + + rasterBand.getBandWidth() = (maxX-origin.getX())*mx + */ + + double mx = rasterBand.getBandWidth()/(maxX-origin.getX()); + double bx = -origin.getX()*mx; + double px = mx*x + bx; + + double my = rasterBand.getBandHeight()/(maxY-origin.getY()); + double by = -origin.getY()*my; + double py = my*y + by; + + /* + if (debug) { + log.debug("p(x, y): " + px + " / " + py); + } + */ + SeRasterConstraint constraint = new SeRasterConstraint(); - constraint.setLevel(maxLevel); - constraint.setBands(bandNumber); - constraint.setEnvelope(pixelX / rasterSize, pixelY / rasterSize, pixelX / rasterSize, pixelY / rasterSize); + constraint.setLevel(0); // best resolution + constraint.setBands(rasterBand.getBandNumber()); + int tx = (int)Math.floor(px / rasterBand.getTileWidth()); + int ty = (int)Math.floor(py / rasterBand.getTileHeight()); + constraint.setEnvelope(tx, ty, tx, ty); + pSeQuery.queryRasterTile(constraint); SeRasterTile tile = row.getRasterTile(); - - if(tile != null){ - log.debug("BAND-ID Tile "+tile.getBandId().longValue()); - log.debug("Pixel "+ tile.getNumPixels()); - log.debug("Column / Row "+tile.getColumnIndex()+" / "+tile.getRowIndex()); - - + + if (tile != null){ double[] tileValues = new double[tile.getNumPixels()]; tileValues = tile.getPixels(tileValues); - lSet.addCol(new ColDefinition("tile",ColDefinition.FLOAT64)); + lSet.addCol(new ColDefinition("tile", ColDefinition.FLOAT64)); Row lBackingRow = new Row(1); - Envelope envelope = new Envelope( - new Coordinate(extent.getMinX(), - extent.getMinY()), - new Coordinate(extent.getMaxX(), - extent.getMaxY())); - RasterObject ro = new RasterObject(envelope, - tileValues, - rasterSize, - tile.getRowIndex(), - tile.getColumnIndex(), - width, - height); - lBackingRow.addObject(ro,0); + + /* + 0 = wx1*mt + bt + tileWidth-EPS = wx2*mt + bt + tileWidth-EPS = mt*(wx2 - wx1) + */ + + double wx1 = (rasterBand.getTileWidth()*tx - bx)/mx; + double wx2 = (rasterBand.getTileWidth()*(tx+1) - bx)/mx; + double mxt = (rasterBand.getTileWidth()-1e-5d)/(wx2 - wx1); + double bxt = -wx1*mxt; + + double wy1 = (rasterBand.getTileHeight()*ty - by)/my; + double wy2 = (rasterBand.getTileHeight()*(ty+1) - by)/my; + double myt = (rasterBand.getTileHeight()-1e-5d)/(wy2 - wy1); + double byt = -wy1*myt; + + RasterObject ro = new RasterObject( + mxt, bxt, + myt, byt, + tile.getColumnIndex(), + tile.getRowIndex(), + tileValues, + rasterBand.getTileWidth(), + rasterBand.getTileHeight()); + lBackingRow.addObject(ro, 0); lSet.addRow(lBackingRow); + + /* + if (debug) { + log.debug("x / y: " + x + " / " + y); + log.debug("wx1: " + wx1); + log.debug("wx1 -> " + (wx1*mxt+bxt)); + log.debug("wx2: " + wx2); + log.debug("wx2 -> " + (wx2*mxt+bxt)); + log.debug("wx2 - wx1: " + Math.abs(wx2-wx1)); + log.debug("wy2 - wy1: " + Math.abs(wy2-wy1)); + log.debug("pix: " + (x*mxt + bxt)+ " " + (y*myt + byt)); + log.debug("requesting tile: " + tx + " / " + ty); + log.debug("got tile: " + tile.getColumnIndex() + " / " + tile.getRowIndex()); + log.debug("tile orig: " + origin.getX() + " / " + origin.getY()); + log.debug("tile width: " + rasterBand.getTileWidth()); + log.debug("tile height: " + rasterBand.getTileHeight()); + log.debug("Rasterdimesion (Pixel) " + + rasterBand.getBandWidth() + " / " + rasterBand.getBandHeight()); + log.debug("BAND-ID "+rasterBand.getId().longValue()); + log.debug("Pixels: "+ tile.getNumPixels()); + log.debug("Column / Row "+tile.getColumnIndex()+" / "+tile.getRowIndex()); + } + */ } - - }else{ - log.debug("The Query doesn't deliver any Information because the Point is located outside the Raster."); + } + else{ + log.debug("The Query doesn't deliver any Information " + + "because the Point is located outside the Raster."); } } catch (Exception e) { log.error(e,e); @@ -613,6 +676,6 @@ public <T> T unwrap(Class<T> iface) throws SQLException { return null; + } - }
--- a/geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java Fri Jan 08 14:40:08 2010 +0000 +++ b/geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java Sat Jan 09 12:41:26 2010 +0000 @@ -8,75 +8,105 @@ import java.io.Serializable; + +import org.apache.log4j.Logger; + /** - * @author Tim Englich <tim.englich@intevation.de> - * + * @author Tim Englich (tim.englich@intevation.de) + * @author Sascha L. Teichmann (sascha.teichmann@intevation.de) */ public class RasterObject implements Serializable { - private int rasterSize; - private int rowIndex; + private static Logger log = Logger.getLogger(RasterObject.class); + + public static final int NEAREST_NEIGHBOR = 0; + public static final int BILINEAR = 1; + + private double mx; + private double bx; + private double my; + private double by; + private int columnIndex; - private int rasterWidth; - private int rasterHeight; - private Envelope envelope = null; - private Envelope rasterEnvelope = null; - - private double[] rasterData = null; - /** - * Constructor - */ - public RasterObject(Envelope envelope , double[] rasterData, - int rasterSize, int rowIndex, - int columnIndex, int rasterWidth, int rasterHeight ){ - this.envelope = envelope; - this.rasterData = rasterData; - this.rasterSize = rasterSize; - this.rowIndex = rowIndex; + private int rowIndex; + + private int tileWidth; + private int tileHeight; + + private double [] rasterData; + + public RasterObject() { + } + + public RasterObject( + double mx, double bx, + double my, double by, + int columnIndex, + int rowIndex, + double [] rasterData, + int tileWidth, + int tileHeight + ){ + this.mx = mx; + this.bx = bx; + this.my = my; + this.by = by; this.columnIndex = columnIndex; - this.rasterHeight = rasterHeight; - this.rasterWidth = rasterWidth; + this.rowIndex = rowIndex; + this.rasterData = rasterData; + this.tileWidth = tileWidth; + this.tileHeight = tileHeight; } - public synchronized Envelope getEnvelope() { - if (this.rasterEnvelope == null){ - double minX = (columnIndex* rasterSize * - (envelope.getWidth() / rasterWidth)+ envelope.getMinX()); - double maxX = minX + (rasterSize * (envelope.getWidth() / rasterWidth)); - - double maxY = envelope.getMaxY() - - (rowIndex * rasterSize * - (envelope.getHeight() / rasterHeight)); - double minY = maxY - - (rasterSize * (envelope.getHeight() / rasterHeight)); - this.rasterEnvelope = new Envelope( - new Coordinate(minX, minY), - new Coordinate(maxX, maxY)); - } - return this.rasterEnvelope; + + public boolean contains(Coordinate coordinate) { + double px = mx*coordinate.x + bx; + double py = my*coordinate.y + by; + return px >= 0d && py >= 0d && px < tileWidth && py < tileHeight; + } + + public final double get(int posX, int posY) { + return rasterData[posY*tileWidth + posX]; + } + + public int getTileWidth() { + return tileWidth; + } + + public int getTileHeight() { + return tileWidth; + } + + public int getColumnIndex() { + return columnIndex; + } + + public int getRowIndex() { + return rowIndex; } - public synchronized double getValue(Coordinate coordinate){ - double dxNature = coordinate.x - envelope.getMinX(); - double dyNature = envelope.getMaxY() - coordinate.y; - double widthNature = envelope.getMaxX()-envelope.getMinX(); - double heightNature = envelope.getMaxY()-envelope.getMinY(); - - int pixelX = (int)Math.round(dxNature * rasterWidth / widthNature); - int pixelY = (int)Math.round(dyNature * rasterHeight/ heightNature); - - - int localPixelX = pixelX - (columnIndex*rasterSize); - int localPixelY = pixelY - (rowIndex*rasterSize); - - int pos = localPixelY * rasterSize + localPixelX; - - - if (pos < rasterData.length){ - return this.rasterData[pos]; - } - - return Double.NaN; + public double getValue(Coordinate coordinate) { + return getValue(coordinate, NEAREST_NEIGHBOR); } + public double getValue(Coordinate coordinate, int interpolationType) { + + double px = mx*coordinate.x + bx; + double py = my*coordinate.y + by; + + if (px < 0d || py < 0d || px >= tileWidth || py >= tileHeight) { + return Double.NaN; + } + + /* + if (interpolationType == BILINEAR) { + // TODO: implement me! + } + */ + + int posX = Math.min(tileWidth-1, (int)Math.round(px)); + int posY = Math.min(tileHeight-1, (int)Math.round(py)); + + return get(posX, posY); + } }