public void writeGrid(GridDatatype grid, Array data, boolean greyScale,
double xStart, double yStart, double xInc,
double yInc, int imageNumber) throws IOException {
int nextStart = 0;
GridCoordSystem gcs = grid.getCoordinateSystem();
// get rid of this when all projections are implemented
if ( !gcs.isLatLon()
&& !(gcs.getProjection() instanceof LambertConformal)
&& !(gcs.getProjection() instanceof Stereographic)
&& !(gcs.getProjection() instanceof Mercator)
// && !(gcs.getProjection() instanceof TransverseMercator)
&& !(gcs.getProjection() instanceof AlbersEqualAreaEllipse)
&& !(gcs.getProjection() instanceof AlbersEqualArea)) {
throw new IllegalArgumentException(
"Must be lat/lon or LambertConformal or Mercator and grid = "
+ gcs.getProjection().getClass().getName());
}
// write the data first
if (greyScale) {
ArrayByte result = replaceMissingValuesAndScale(grid, data);
nextStart = geotiff.writeData((byte[]) result.getStorage(),
imageNumber);
} else {
ArrayFloat result = replaceMissingValues(grid, data);
nextStart = geotiff.writeData((float[]) result.getStorage(),
imageNumber);
}
// set the width and the height
int elemSize = greyScale
? 1
: 4;
int height = data.getShape()[0]; // Y
int width = data.getShape()[1]; // X
int size = elemSize * height * width; // size in bytes
geotiff.addTag(new IFDEntry(Tag.ImageWidth,
FieldType.SHORT).setValue(width));
geotiff.addTag(new IFDEntry(Tag.ImageLength,
FieldType.SHORT).setValue(height));
// set the multiple images tag
int ff = 1 << 1;
int page = imageNumber - 1;
geotiff.addTag(new IFDEntry(Tag.NewSubfileType,
FieldType.SHORT).setValue(ff));
geotiff.addTag(new IFDEntry(Tag.PageNumber,
FieldType.SHORT).setValue(page, 2));
// just make it all one big "row"
geotiff.addTag(new IFDEntry(Tag.RowsPerStrip,
FieldType.SHORT).setValue(1)); //height));
// the following changes to make it viewable in ARCMAP
/*
geotiff.addTag( new IFDEntry(Tag.StripByteCounts, FieldType.LONG).setValue( size));
// data starts here, header is written at the end
if( imageNumber == 1 )
geotiff.addTag( new IFDEntry(Tag.StripOffsets, FieldType.LONG).setValue( 8));
else
geotiff.addTag( new IFDEntry(Tag.StripOffsets, FieldType.LONG).setValue(nextStart));
*/
int[] soffset = new int[height];
int[] sbytecount = new int[height];
if (imageNumber == 1) {
soffset[0] = 8;
} else {
soffset[0] = nextStart;
}
sbytecount[0] = width * elemSize;
for (int i = 1; i < height; i++) {
soffset[i] = soffset[i - 1] + width * elemSize;
sbytecount[i] = width * elemSize;
}
geotiff.addTag(new IFDEntry(Tag.StripByteCounts, FieldType.LONG,
width).setValue(sbytecount));
geotiff.addTag(new IFDEntry(Tag.StripOffsets, FieldType.LONG,
width).setValue(soffset));
// standard tags
geotiff.addTag(new IFDEntry(Tag.Orientation,
FieldType.SHORT).setValue(1));
geotiff.addTag(new IFDEntry(Tag.Compression,
FieldType.SHORT).setValue(1)); // no compression
geotiff.addTag(new IFDEntry(Tag.Software,
FieldType.ASCII).setValue("nc2geotiff"));
geotiff.addTag(new IFDEntry(Tag.PhotometricInterpretation,
FieldType.SHORT).setValue(1)); // black is zero : not used?
geotiff.addTag(new IFDEntry(Tag.PlanarConfiguration,
FieldType.SHORT).setValue(1));
if (greyScale) {
// standard tags for Greyscale images ( see TIFF spec, section 4)
geotiff.addTag(new IFDEntry(Tag.BitsPerSample,
FieldType.SHORT).setValue(8)); // 8 bits per sample
geotiff.addTag(new IFDEntry(Tag.SamplesPerPixel,
FieldType.SHORT).setValue(1));
geotiff.addTag(new IFDEntry(Tag.XResolution,
FieldType.RATIONAL).setValue(1, 1));
geotiff.addTag(new IFDEntry(Tag.YResolution,
FieldType.RATIONAL).setValue(1, 1));
geotiff.addTag(new IFDEntry(Tag.ResolutionUnit,
FieldType.SHORT).setValue(1));
} else {
// standard tags for SampleFormat ( see TIFF spec, section 19)
geotiff.addTag(new IFDEntry(Tag.BitsPerSample,
FieldType.SHORT).setValue(32)); // 32 bits per sample
geotiff.addTag(new IFDEntry(Tag.SampleFormat,
FieldType.SHORT).setValue(3)); // Sample Format
geotiff.addTag(new IFDEntry(Tag.SamplesPerPixel,
FieldType.SHORT).setValue(1));
MAMath.MinMax dataMinMax = grid.getMinMaxSkipMissingData(data);
float min = (float) (dataMinMax.min);
float max = (float) (dataMinMax.max);
geotiff.addTag(new IFDEntry(Tag.SMinSampleValue,
FieldType.FLOAT).setValue(min));
geotiff.addTag(new IFDEntry(Tag.SMaxSampleValue,
FieldType.FLOAT).setValue(max));
geotiff.addTag(new IFDEntry(Tag.GDALNoData,
FieldType.FLOAT).setValue(min - 1.f));
}
/*
geotiff.addTag( new IFDEntry(Tag.Geo_ModelPixelScale, FieldType.DOUBLE).setValue(
new double[] {5.0, 2.5, 0.0} ));
geotiff.addTag( new IFDEntry(Tag.Geo_ModelTiepoint, FieldType.DOUBLE).setValue(
new double[] {0.0, 0.0, 0.0, -180.0, 90.0, 0.0 } ));
// new double[] {0.0, 0.0, 0.0, 183.0, 90.0, 0.0} ));
IFDEntry ifd = new IFDEntry(Tag.Geo_KeyDirectory, FieldType.SHORT).setValue(
new int[] {1, 1, 0, 4, 1024, 0, 1, 2, 1025, 0, 1, 1, 2048, 0, 1, 4326, 2054, 0, 1, 9102} );
geotiff.addTag( ifd);
*/
// set the transformation from projection to pixel, add tie point tag
geotiff.setTransform(xStart, yStart, xInc, yInc);
if (gcs.isLatLon()) {
addLatLonTags();
} else if (gcs.getProjection() instanceof LambertConformal) {
addLambertConformalTags((LambertConformal) gcs.getProjection(),
xStart, yStart);
} else if (gcs.getProjection() instanceof Stereographic) {
addPolarStereographicTags((Stereographic) gcs.getProjection(),
xStart, yStart);
} else if (gcs.getProjection() instanceof Mercator) {
addMercatorTags((Mercator) gcs.getProjection());
} else if (gcs.getProjection() instanceof TransverseMercator) {
addTransverseMercatorTags(
(TransverseMercator) gcs.getProjection());
} else if (gcs.getProjection() instanceof AlbersEqualArea) {
addAlbersEqualAreaTags((AlbersEqualArea) gcs.getProjection());
} else if (gcs.getProjection() instanceof AlbersEqualAreaEllipse) {
addAlbersEqualAreaEllipseTags((AlbersEqualAreaEllipse) gcs.getProjection());
}
geotiff.writeMetadata(imageNumber);
//geotiff.close();