Hello all,
I received a question about the rasterization code (Features to Grid
Coverage) and upgraded the code to Geotools 2.5.2. I've attached the
main class (FeatureRasterizer.java) and a test class
(TestFeatureRasterizer.java).
Here are the geotools jars needed to run this:
geoapi-2.2-M1.jar
gt-api-2.5.2.jar
gt-coverage-2.5.2.jar
gt-epsg-wkt-2.5.2.jar
gt-main-2.5.2.jar
gt-metadata-2.5.2.jar
gt-referencing-2.5.2.jar
gt-render-2.5.2.jar
gt-shapefile-2.5.2.jar
jai_codec-1.1.3.jar
jai_imageio-1.1.jar
jsr-275-1.0-beta-2.jar
jts-1.9.jar
vecmath-1.3.1.jar
Steve
Steve Ansari said the following on 4/1/2008 4:46 PM:
> Thanks Ferdinando! It looks great.
>
> Steve
>
>
> Ferdinando Villa said the following on 4/1/2008 4:18 PM:
>> Actually Steve's code works beautifully and does not seem to require much
>> modernization! I haven't tested it with a subsampled envelope yet, but it
>> spits out 4000x2000 rasters from big shapefiles in seconds. I only commented
>> out some (not all) debug printouts and added a function to produce a
>> GridCoverage2D directly:
>>
>> public GridCoverage2D rasterize(String name, FeatureCollection fc,
>> String attributeName, ReferencedEnvelope env)
>>
>> you can pass null as the envelope and it will use the boundaries of the
>> feature collection. I also internalized the FeatureRasterizerException so
>> the code is now self-contained as long as geotools is in the classpath. It
>> does not use any deprecated API - that's the extent of my modernization
>> capabilities :)
>>
>> If I find bugs along the way I'll try to fix them - progress can be followed
>> at
>> http://ecoinformatics.uvm.edu/crucible/browse/ThinklabGeospacePlugin/org/int
>> egratedmodelling/geospace/gis/FeatureRasterizer.java.
>>
>> Bottom line: Thanks Steve!
>>
>> Cheers,
>> ferdinando
>>
>>
>> -----Original Message-----
>> From: Jody Garnett [mailto:[email protected]]
>> Sent: Tuesday, April 01, 2008 3:28 PM
>> To: Ferdinando Villa
>> Cc: 'Steve Ansari'; [email protected]
>> Subject: Re: [Geotools-gt2-users] Vector to Raster conversion
>>
>> Hey guys; we are setting up a "Process" module in geotools (see this
>> weeks IRC). After you do your "testing and modernize" do you want to
>> send it our way? It woudl be great to have examples of different kinds
>> of processes to make sure the API works out okay.
>> Jody
>>
>> ------------------------------------------------------------------------
>>
>> -------------------------------------------------------------------------
>> Check out the new SourceForge.net Marketplace.
>> It's the best place to buy or sell services for
>> just about anything Open Source.
>> http://ad.doubleclick.net/clk;164216239;13503038;w?http://sf.net/marketplace
>> ------------------------------------------------------------------------
>>
>> _______________________________________________
>> Geotools-gt2-users mailing list
>> [email protected]
>> https://lists.sourceforge.net/lists/listinfo/geotools-gt2-users
>>
>
> --
> Steve Ansari
> Physical Scientist
> NOAA's National Climatic Data Center
> Veach-Baley Federal Building
> 151 Patton Avenue
> Asheville, NC 28801
> Ph: 828-271-4611
> Fax: 828-271-4022
--
Steve Ansari
Physical Scientist
NOAA's National Climatic Data Center
Veach-Baley Federal Building
151 Patton Avenue
Asheville, NC 28801
Ph: 828-271-4611
Fax: 828-271-4328
package gov.noaa.ncdc.geotools;
/**
* NOAA's National Climatic Data Center
* NOAA/NESDIS/NCDC
* 151 Patton Ave, Asheville, NC 28801
*
* THIS SOFTWARE AND ITS DOCUMENTATION ARE CONSIDERED TO BE IN THE
* PUBLIC DOMAIN AND THUS ARE AVAILABLE FOR UNRESTRICTED PUBLIC USE.
* THEY ARE FURNISHED "AS IS." THE AUTHORS, THE UNITED STATES GOVERNMENT, ITS
* INSTRUMENTALITIES, OFFICERS, EMPLOYEES, AND AGENTS MAKE NO WARRANTY,
* EXPRESS OR IMPLIED, AS TO THE USEFULNESS OF THE SOFTWARE AND
* DOCUMENTATION FOR ANY PURPOSE. THEY ASSUME NO RESPONSIBILITY (1)
* FOR THE USE OF THE SOFTWARE AND DOCUMENTATION; OR (2) TO PROVIDE
* TECHNICAL SUPPORT TO USERS.
*/
import java.awt.AlphaComposite;
import java.awt.Color;
import java.awt.Graphics2D;
import java.awt.GraphicsEnvironment;
import java.awt.image.BufferedImage;
import java.awt.image.DataBuffer;
import java.awt.image.WritableRaster;
import java.nio.ByteBuffer;
import javax.media.jai.RasterFactory;
import org.geotools.coverage.grid.GridCoverage2D;
import org.geotools.coverage.grid.GridCoverageFactory;
import org.geotools.feature.FeatureCollection;
import org.geotools.feature.FeatureIterator;
import org.geotools.geometry.jts.ReferencedEnvelope;
import org.opengis.feature.Feature;
import com.vividsolutions.jts.geom.Coordinate;
import com.vividsolutions.jts.geom.Envelope;
import com.vividsolutions.jts.geom.Geometry;
import com.vividsolutions.jts.geom.GeometryFactory;
import com.vividsolutions.jts.geom.LineString;
import com.vividsolutions.jts.geom.LinearRing;
import com.vividsolutions.jts.geom.MultiLineString;
import com.vividsolutions.jts.geom.MultiPoint;
import com.vividsolutions.jts.geom.MultiPolygon;
import com.vividsolutions.jts.geom.Point;
import com.vividsolutions.jts.geom.Polygon;
/**
* Rasterize features onto a WritableRaster object using Java 2D
Graphics/BufferedImage.
*
* @author steve.ansari
* @author Ferdinando Villa
* @created March 20, 2008
*/
public class FeatureRasterizer {
public class FeatureRasterizerException extends Exception {
/**
* Constructor with message argument.
*
* @param message Reason for the exception being thrown
*/
public FeatureRasterizerException(String message) {
super(message);
}
}
private int height;
private int width;
private double noDataValue;
private WritableRaster raster = null;
private BufferedImage bimage = null;
private Graphics2D graphics = null;
private java.awt.geom.Rectangle2D.Double bounds;
private double cellsize;
private double minAttValue = 999999999;
private double maxAttValue = -999999999;
// Declare these as global
private int[] coordGridX = new int[3500];
private int[] coordGridY = new int[3500];
private float value;
private boolean emptyGrid = false;
private Geometry extentGeometry;
private GeometryFactory geoFactory = new GeometryFactory();
private String attributeName = "value";
public static GridCoverageFactory rasterFactory = new GridCoverageFactory();
private double xInterval;
private double yInterval;
// Any change in height, width or no_data values will cause
// the raster to 'reset' at the next call to .rasterize(...)
private boolean resetRaster = false;
/**
*Constructor for the FeatureRasterizer object
*
* @exception FeatureRasterizerException Description of the Exception
*/
public FeatureRasterizer() {
this(800, 800, -999.0f);
}
/**
* Constructor for the FeatureRasterizer object - will use default 800x800
raster
*
* @param noData No Data value for raster
* @exception FeatureRasterizerException Description of the Exception
*/
public FeatureRasterizer(float noData) {
this(800, 800, noData);
}
/**
* Constructor for the FeatureRasterizer object. No Data value defaults to
-999.0
*
* @param height Height of raster (number of grid
cells)
* @param width Width of raster (number of grid
cells)
*/
public FeatureRasterizer(int height, int width) {
this(height, width, -999.0f);
}
/**
* Constructor for the FeatureRasterizer object
*
* @param height Height of raster (number of grid
cells)
* @param width Width of raster (number of grid
cells)
* @param noData No Data value for raster
*/
public FeatureRasterizer(int height, int width, float noData) {
this.height = height;
this.width = width;
this.noDataValue = noData;
raster = RasterFactory.createBandedRaster(DataBuffer.TYPE_FLOAT,
width, height, 1, null);
bimage = new BufferedImage(width, height, BufferedImage.TYPE_INT_ARGB);
bimage.setAccelerationPriority(1.0f);
GraphicsEnvironment ge =
GraphicsEnvironment.getLocalGraphicsEnvironment();
// System.out.println("IMAGE ACCELERATED?
"+bimage.getCapabilities(ge.getDefaultScreenDevice().getDefaultConfiguration()).isAccelerated());
graphics = bimage.createGraphics();
graphics.setPaintMode();
graphics.setComposite(AlphaComposite.Src);
}
public GridCoverage2D rasterize(String name, FeatureCollection fc, String
attributeName, ReferencedEnvelope env) throws FeatureRasterizerException {
if (raster == null) {
WritableRaster raster =
RasterFactory.createBandedRaster(
DataBuffer.TYPE_FLOAT,
this.width,
this.height,
1,
null);
setWritableRaster(raster);
}
clearRaster();
if (env == null) {
rasterize(fc, attributeName);
/*
* Use full envelope from feature collection
* TODO check if we need to use a buffer like in Steve's code above
*/
env = fc.getBounds();
} else {
/*
* TODO check if we need to use a buffer like in Steve's code above
*/
java.awt.geom.Rectangle2D.Double box =
new java.awt.geom.Rectangle2D.Double(
env.getMinX(),
env.getMinY(),
env.getWidth(),
env.getHeight());
rasterize(fc, box, attributeName);
}
GridCoverage2D coverage =
rasterFactory.create(name, raster, env);
return coverage;
}
/**
* Gets the raster attribute of the FeatureRasterizer object
* Processes data from the FeatureCollection and approximates onto a
Raster Grid.
*
* @param fc Feature Collection with features
to rasterize.
* @param attributeName Name of attribute from feature
collection to provide as the cell value.
* @exception FeatureRasterizerException An error when rasterizing the
data
*/
public void rasterize(FeatureCollection fc, String attributeName)
throws FeatureRasterizerException {
// calculate variable resolution bounds that fit around feature
collection
double edgeBuffer = 0.001;
double x = fc.getBounds().getMinX() - edgeBuffer;
double y = fc.getBounds().getMinY() - edgeBuffer;
double width = fc.getBounds().getWidth() + edgeBuffer * 2;
double height = fc.getBounds().getHeight() + edgeBuffer * 2;
java.awt.geom.Rectangle2D.Double bounds = new
java.awt.geom.Rectangle2D.Double(x, y, width, height);
System.out.println("BOUNDS: "+bounds);
System.out.println("FCBNDS: "+fc.getBounds());
rasterize(fc, bounds, attributeName);
}
/**
* Gets the raster attribute of the FeatureRasterizer object
* Processes data from the FeatureCollection and approximates onto a
Raster Grid.
*
* @param fc Description of the Parameter
* @param bounds Description of the Parameter
* @param attributeName Name of attribute from feature
collection to provide as the cell value.
* @exception FeatureRasterizerException An error when rasterizing the
data
*/
public void rasterize(FeatureCollection fc,
java.awt.geom.Rectangle2D.Double bounds, String attributeName)
throws FeatureRasterizerException {
this.attributeName = attributeName;
// Check if we need to change the underlying raster
if (resetRaster) {
raster = RasterFactory.createBandedRaster(DataBuffer.TYPE_FLOAT,
width, height, 1, null);
bimage = new BufferedImage(width, height,
BufferedImage.TYPE_INT_ARGB);
bimage.setAccelerationPriority(1.0f);
GraphicsEnvironment ge =
GraphicsEnvironment.getLocalGraphicsEnvironment();
// System.out.println("IMAGE ACCELERATED?
"+bimage.getCapabilities(ge.getDefaultScreenDevice().getDefaultConfiguration()).isAccelerated());
graphics = bimage.createGraphics();
graphics.setPaintMode();
graphics.setComposite(AlphaComposite.Src);
resetRaster = false;
//System.out.println("---------------- RESETING FeatureRasterizer
WritableRaster OBJECT -------------------- ");
}
// initialize raster to NoData value
clearRaster();
setBounds(bounds);
// TODO - change method calls to account for a switch to control if
rasterizer should work if vis bounds > feature bounds
// All the data should start in the lower left corner. Don't export
what we don't need.
double ratio = bounds.height / bounds.width;
int ncols;
int nrows;
if (ratio < 1) {
// wider than tall
nrows = (int) (ratio * height);
ncols = width;
}
else {
nrows = height;
ncols = (int) (height / ratio);
}
//System.out.println("1 --- WIDTH: " + ncols + " HEIGHT: " + nrows);
FeatureIterator fci = fc.features();
Feature feature;
while (fci.hasNext()) {
feature = fci.next();
addFeature(feature);
}
close();
}
/**
* Implementation the StreamingProcess interface. Rasterize a single
feature and
* update current WriteableRaster using the current settings.
*
* @param feature The feature to rasterize and add to current
WritableRaster
*/
public void addFeature(Feature feature) {
//System.out.println("rasterizer - processing feature: "+ attributeName
+ " -- " + feature);
try {
value =
Float.parseFloat(feature.getProperty(attributeName).getValue().toString());
if (value > maxAttValue) { maxAttValue = value; }
if (value < minAttValue) { minAttValue = value; }
} catch (Exception e) {
e.printStackTrace();
System.err.println("THE FEATURE COULD NOT BE RASTERIZED BASED ON
THE '"+attributeName+
"' ATTRIBUTE VALUE OF
'"+feature.getProperty(attributeName).getValue().toString()+"'");
return;
}
int rgbVal = floatBitsToInt(value);
graphics.setColor(new Color(rgbVal, true));
// Extract polygon and rasterize!
Geometry geometry =
(Geometry)feature.getDefaultGeometryProperty().getValue();
if (geometry.intersects(extentGeometry)) {
if (geometry.getClass().equals(MultiPolygon.class)) {
MultiPolygon mp = (MultiPolygon)geometry;
for (int n=0; n<mp.getNumGeometries(); n++) {
drawGeometry(mp.getGeometryN(n));
}
}
else if (geometry.getClass().equals(MultiLineString.class)) {
MultiLineString mp = (MultiLineString)geometry;
for (int n=0; n<mp.getNumGeometries(); n++) {
drawGeometry(mp.getGeometryN(n));
}
}
else if (geometry.getClass().equals(MultiPoint.class)) {
MultiPoint mp = (MultiPoint)geometry;
for (int n=0; n<mp.getNumGeometries(); n++) {
drawGeometry(mp.getGeometryN(n));
}
}
else {
drawGeometry(geometry);
}
}
}
/**
* Implementation the StreamingProcess interface - this copies values from
BufferedImage RGB to WritableRaster of floats.
*/
public void close() {
for (int i = 0; i < width; i++) {
for (int j = 0; j < height; j++) {
double val = Float.intBitsToFloat(bimage.getRGB(i, j));
raster.setSample(i, j, 0, val);
}
}
}
private void drawGeometry(Geometry geometry) {
Coordinate[] coords = geometry.getCoordinates();
// enlarge if needed
if (coords.length > coordGridX.length) {
coordGridX = new int[coords.length];
coordGridY = new int[coords.length];
}
// Clear Array
for (int i = 0; i < coords.length; i++) {
coordGridX[i] = -1;
}
for (int i = 0; i < coords.length; i++) {
coordGridY[i] = -1;
}
// Go through coordinate array in order received (clockwise)
for (int n = 0; n < coords.length; n++) {
coordGridX[n] = (int) (((coords[n].x - bounds.x) / xInterval));
coordGridY[n] = (int) (((coords[n].y - bounds.y) / yInterval));
coordGridY[n] = bimage.getHeight() - coordGridY[n];
}
if (geometry.getClass().equals(Polygon.class)) {
graphics.fillPolygon(coordGridX, coordGridY, coords.length);
}
else if (geometry.getClass().equals(LinearRing.class)) {
graphics.drawPolyline(coordGridX, coordGridY, coords.length);
}
else if (geometry.getClass().equals(LineString.class)) {
graphics.drawPolyline(coordGridX, coordGridY, coords.length);
}
else if (geometry.getClass().equals(Point.class)) {
graphics.drawPolyline(coordGridX, coordGridY, coords.length);
}
}
/**
* Gets the emptyGrid attribute of the FeatureRasterizer object
*
* @return The emptyGrid value
*/
public boolean isEmptyGrid() {
return emptyGrid;
}
/**
* Gets the writableRaster attribute of the FeatureRasterizer object
*
* @return The writableRaster value
*/
public WritableRaster getWritableRaster() {
return raster;
}
/**
* Sets the writableRaster attribute of the FeatureRasterizer object
*
* @param raster The new writableRaster value
*/
public void setWritableRaster(WritableRaster raster) {
this.raster = raster;
}
/**
* Gets the bounds attribute of the FeatureRasterizer object
*
* @return The bounds value
*/
public java.awt.geom.Rectangle2D.Double getBounds() {
return bounds;
}
/**
* Sets the bounds for the Rasterizer
*
* @return The bounds value
*/
public void setBounds(java.awt.geom.Rectangle2D.Double bounds) {
this.bounds = bounds;
xInterval = bounds.width / (double) width;
yInterval = bounds.height / (double) height;
System.out.println("xInterval: " + xInterval + " yInterval: " +
yInterval);
if (xInterval > yInterval) {
yInterval = xInterval;
}
if (yInterval > xInterval) {
xInterval = yInterval;
}
cellsize = yInterval;
// Clip geometries to the provided bounds
// Create extent geometry
Envelope env = new Envelope(
bounds.getX(),
bounds.getX() + bounds.getWidth(),
bounds.getY(),
bounds.getY() + bounds.getHeight()
);
extentGeometry = geoFactory.toGeometry(env);
}
/**
* Sets the entire raster to NoData
*/
public void clearRaster() {
// System.out.println("CLEARING RASTER");
minAttValue = 999999999;
maxAttValue = -999999999;
// initialize raster to NoData value
for (int i = 0; i < width; i++) {
for (int j = 0; j < height; j++) {
raster.setSample(i, j, 0, noDataValue);
bimage.setRGB(i, j, floatBitsToInt((float)noDataValue));
}
}
}
/**
* Get the current attribute to use as the grid cell values.
*/
public String getAttName() {
return attributeName;
}
/**
* Sets the current attribute to use as the grid cell values.
*/
public void setAttName(String attName) {
this.attributeName = attName;
}
/**
* Gets the cellsize attribute of the FeatureRasterizer object
*
* @return The cellsize value
*/
public double getCellsize() {
return cellsize;
}
public double getNoDataValue() {
return noDataValue;
}
public void setNoDataValue(double noData) {
if (noData != noDataValue) {
resetRaster = true;
}
this.noDataValue = noData;
}
public int getHeight() {
return height;
}
public void setHeight(int height) {
if (height != height) {
resetRaster = true;
}
this.height = height;
}
public int getWidth() {
return width;
}
public void setWidth(int width) {
if (width != width) {
resetRaster = true;
}
this.width = width;
}
public double getMinAttValue() {
return minAttValue;
}
public double getMaxAttValue() {
return maxAttValue;
}
private static int floatBitsToInt(float f) {
ByteBuffer conv = ByteBuffer.allocate(4);
conv.putFloat(0, f);
return conv.getInt(0);
}
public String toString() {
return "FEATURE RASTERIZER: WIDTH="+width+" , HEIGHT="+height+" ,
NODATA="+noDataValue;
}
}
package gov.noaa.ncdc.geotools;
import gov.noaa.ncdc.geotools.FeatureRasterizerSteve.FeatureRasterizerException;
import java.awt.geom.Rectangle2D;
import java.io.File;
import java.io.IOException;
import java.net.URL;
import java.util.NoSuchElementException;
import javax.swing.JFrame;
import org.geotools.coverage.grid.GridCoverage2D;
import org.geotools.coverage.grid.GridCoverageFactory;
import org.geotools.data.FeatureSource;
import org.geotools.data.shapefile.ShapefileDataStore;
import org.geotools.feature.FeatureCollection;
import org.geotools.feature.FeatureIterator;
import org.geotools.geometry.jts.ReferencedEnvelope;
import org.geotools.referencing.CRS;
import org.opengis.geometry.MismatchedDimensionException;
import org.opengis.referencing.FactoryException;
import org.opengis.referencing.NoSuchAuthorityCodeException;
import com.sun.media.jai.widget.DisplayJAI;
public class TestFeatureRasterizer {
public static void main(String[] args) {
try {
testRasterizer();
} catch (Exception e) {
e.printStackTrace();
}
}
public static void testRasterizer() throws IOException,
NoSuchElementException, FeatureRasterizerException,
MismatchedDimensionException, NoSuchAuthorityCodeException, FactoryException {
URL url = new File("H:\\ESRI\\usa\\counties.shp").toURI().toURL();
ShapefileDataStore ds = new ShapefileDataStore(url);
FeatureSource fs = ds.getFeatureSource("counties");
FeatureRasterizerSteve rasterizer = new FeatureRasterizerSteve(800,
800, -999.0f);
rasterizer.setAttName("POP2000");
// Envelope env = fs.getBounds();
// Rectangle2D.Double bounds = new Rectangle2D.Double(env.getMinX(),
env.getMinY(), env.getWidth(), env.getHeight());
Rectangle2D.Double bounds = new Rectangle2D.Double(-120.0, 20.0, 40.0,
20.0);
rasterizer.setBounds(bounds);
FeatureIterator fi = fs.getFeatures().features();
while (fi.hasNext()) {
rasterizer.addFeature(fi.next());
}
rasterizer.close();
GridCoverageFactory gcFactory = new GridCoverageFactory();
GridCoverage2D gc = gcFactory.create("TEST1",
rasterizer.getWritableRaster(),
new ReferencedEnvelope(bounds.getMinX(), bounds.getMaxX(),
bounds.getMinY(), bounds.getMaxY(),
CRS.decode("EPSG:4326")));
DisplayJAI display = new DisplayJAI(gc.getRenderedImage());
JFrame frame = new JFrame();
frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
frame.getContentPane().add(display);
frame.pack();
frame.setVisible(true);
FeatureCollection fc = fs.getFeatures();
System.out.println("POP2001 feature collection size: "+fc.size());
rasterizer.rasterize(fc, "POP2001");
GridCoverage2D gc2 = gcFactory.create("TEST1",
rasterizer.getWritableRaster(),
new ReferencedEnvelope(bounds.getMinX(), bounds.getMaxX(),
bounds.getMinY(), bounds.getMaxY(),
CRS.decode("EPSG:4326")));
DisplayJAI display2 = new DisplayJAI(gc2.getRenderedImage());
JFrame frame2 = new JFrame();
frame2.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
frame2.getContentPane().add(display2);
frame2.pack();
frame2.setVisible(true);
}
}
------------------------------------------------------------------------------
This SF.net email is sponsored by:
SourcForge Community
SourceForge wants to tell your story.
http://p.sf.net/sfu/sf-spreadtheword
_______________________________________________
Geotools-gt2-users mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/geotools-gt2-users
