RE: [Fwd: multi band tiff from netCDF] - done -ish

Hi,

I have added geotiff tags, I have hacked this code into a class that you can 
use (it has a main method)
but it isn't particulary nicely written (I have spent two days going through 
the JAI source code to
understand it!!)

I have currently restricted the data to WGS84(LL) and floating point variables, 
but this is a restriction on my code, not the method, so it would be fairly 
easy to extend this.  I have commented
out the section that added the time stamp from the netCDF file since it will 
vary from netCDF to netCDF
but it is there to look at.

class usage is java TestGeotiffWriter netCDFFile varName 
and the result is the netCDFFile name with the extension changed to .tif.


Hope this helps,

Norman Barker

p.s. (apologies top posted again, since not sure whether this newsgroup accepts 
attachments and (long) file below)


import java.awt.image.ColorModel;
import java.awt.image.DataBuffer;
import java.awt.image.DataBufferFloat;
import java.awt.image.Raster;
import java.awt.image.SampleModel;
import java.io.File;
import java.io.FileOutputStream;
import java.util.Iterator;
import java.util.List;

import javax.media.jai.JAI;
import javax.media.jai.ParameterBlockJAI;
import javax.media.jai.PlanarImage;
import javax.media.jai.TiledImage;

import ucar.ma2.Array;
import ucar.nc2.DataType;
import ucar.nc2.dataset.CoordinateAxis1D;
import ucar.nc2.dataset.NetcdfDataset;
import ucar.nc2.dataset.grid.GeoGrid;
import ucar.nc2.dataset.grid.GridCoordSys;
import ucar.nc2.dataset.grid.GridDataset;
import ucar.unidata.geoloc.ProjectionPoint;
import ucar.unidata.geoloc.ProjectionRect;

import com.rsinc.data.IngestFile;
import com.sun.media.imageio.plugins.tiff.BaselineTIFFTagSet;
import com.sun.media.imageio.plugins.tiff.GeoTIFFTagSet;
import com.sun.media.jai.codec.ImageCodec;
import com.sun.media.jai.codec.ImageEncoder;
import com.sun.media.jai.codec.TIFFEncodeParam;
import com.sun.media.jai.codec.TIFFField;
import com.sun.media.jai.codecimpl.util.RasterFactory;

public class TestGeotiffWriter {
        public static String EXTENSION = "tif";

        public void convertToGeotiff(String ncdfFilename, String varName) {

                String outGTif = String.copyValueOf(ncdfFilename.toCharArray(), 
0,
                                ncdfFilename.lastIndexOf('.') + 1);
                outGTif = outGTif.concat(TestGeotiffWriter.EXTENSION);
                System.out.println("converting to " + outGTif);

                try {
                        NetcdfDataset ncDataset = 
NetcdfDataset.openDataset(ncdfFilename);

                        GridDataset gSet = new GridDataset(ncDataset);
                        GeoGrid grid = 
gSet.findGridByName(IngestFile.GRID_NAME);
                        GridCoordSys gcs = grid.getCoordinateSystem();

                        if (!gcs.isLatLon()) {
                                System.out.println("Only converting WGS84 (LL) 
at the mo");
                                return;
                        }

                        ProjectionRect bb = gcs.getBoundingBox();
                        ProjectionPoint minPt = bb.getMinPoint();
                        ProjectionPoint maxPt = bb.getMaxPoint();

                        // find the actual tie points of the data
                        double tlX = minPt.getX();
                        double tlY = maxPt.getY();

                        // lat and lon axes seem to get swapped so try this hack
                        // to find res in x and y
                        CoordinateAxis1D xaxis = (CoordinateAxis1D) 
gcs.getXHorizAxis();
                        CoordinateAxis1D yaxis = (CoordinateAxis1D) 
gcs.getYHorizAxis();
                        double xStart = xaxis.getCoordValue(0);
                        double yStart = yaxis.getCoordValue(0);
                        double xEnd = xaxis
                                        
.getCoordValue(grid.getXDimension().getLength() - 1);
                        double yEnd = yaxis
                                        
.getCoordValue(grid.getYDimension().getLength() - 1);
                        double xScale = (xEnd - xStart) / xaxis.getSize();
                        double yScale = (yEnd - yStart) / yaxis.getSize();

                        List list = ncDataset.getCoordinateAxes();
                        Iterator itr = list.iterator();
            boolean found = false;
                        while (itr.hasNext()) {
                                CoordinateAxis1D dim = (CoordinateAxis1D) 
itr.next();
                                String dimName = dim.getName();

                                if (dimName.equals(varName)) {
                                        found = true;
                                        int nBands = (int) dim.getSize();

                                        // do some checks before progressing
                                        if (dim.isUnlimited()) {
                                                // we can't do anything, log
                                                System.out
                                                                
.println("unlimited dimension can write infinite no of bands.");
                                                return;
                                        }

                                        if (dim.getDataType() != 
DataType.FLOAT) {
                                                System.out
                                                                .println("only 
converting variables of floating point at the moment");
                                                return;
                                        }

                                        int width = 
grid.getXDimension().getLength();
                                        int height = 
grid.getYDimension().getLength();

                                        ParameterBlockJAI pbjai = new 
ParameterBlockJAI("bandmerge");

                                        for (int i = 0; i < nBands; i++) {
                                                Array arr = 
grid.readDataSlice(-1, i, -1, -1);
                                                float[] data = (float[]) 
arr.getStorage();
                                                DataBufferFloat dbuffer = new 
DataBufferFloat(data,
                                                                width * height);

                                                // create a float data sample 
model
                                                SampleModel sampleModel = 
RasterFactory
                                                                
.createBandedSampleModel(DataBuffer.TYPE_FLOAT,
                                                                                
width, height, 1);

                                                // create a compatible colour 
model
                                                ColorModel colorModel = 
PlanarImage
                                                                
.createColorModel(sampleModel);

                                                // create a writable raster
                                                Raster raster = 
RasterFactory.createWritableRaster(
                                                                sampleModel, 
dbuffer, null);

                                                // create a tiled image using 
the float sample model
                                                TiledImage tiledImage = new 
TiledImage(0, 0, width,
                                                                height, 0, 0, 
sampleModel, colorModel);

                                                // set the data of the tiled 
image to be the raster
                                                tiledImage.setData(raster);

                                                pbjai.setSource(tiledImage, i);
                                        }

                                        // write to file
                                        PlanarImage result = 
JAI.create("bandmerge", pbjai, null);

                                        File file = new File(outGTif);
                                        FileOutputStream out = new 
FileOutputStream(file);

                                        //Then we create the image encoder and 
encode the buffered image in TIFF.

                                        // write geotiff information, in our 
case WGS84 (LL)
                                        TIFFEncodeParam params = new 
TIFFEncodeParam();

                                        // only dealing with Lat Lon (WGS84)
                                        if (gcs.isLatLon()) {
                                                String[] desc = { 
IngestFile.DIM_NAME + " Bands" };
                                                TIFFField descField = new 
TIFFField(
                                                                
BaselineTIFFTagSet.TAG_IMAGE_DESCRIPTION,
                                                                
TIFFField.TIFF_ASCII, 1, desc);

                                                double[] tiePts = { 0, 0, 0, 
tlX, tlY, 0 };

                                                TIFFField tiePtField = new 
TIFFField(
                                                                
GeoTIFFTagSet.TAG_MODEL_TIE_POINT,
                                                                
TIFFField.TIFF_DOUBLE, tiePts.length, tiePts);

                                                double[] pixelScales = { 
xScale, yScale, 0 };
                                                TIFFField pixelSclField = new 
TIFFField(
                                                                
GeoTIFFTagSet.TAG_MODEL_PIXEL_SCALE,
                                                                
TIFFField.TIFF_DOUBLE, pixelScales.length,
                                                                pixelScales);

                                                char[] geoKeyArr = { 1, 1, 2, 
3, // version, key major, key minor, no. keys
                                                                1024, 0, 1, 2, 
// GTModelTypeGeoKey 1024
                                                                1025, 0, 1, 1, 
// GTRasterTypeGeoKey 1025
                                                                2048, 0, 1, 
4326 // GeographicTypeGeoKey 2048
                                                };
                                                TIFFField geoKeyDir = new 
TIFFField(
                                                                
GeoTIFFTagSet.TAG_GEO_KEY_DIRECTORY,
                                                                
TIFFField.TIFF_SHORT, geoKeyArr.length,
                                                                (Object) 
geoKeyArr);

                                                TIFFField dateField = null;
                                                // get the first time stamp in 
the file
                                                /*                              
                        if (grid.getTimes() != null)
                                                 {
                                                 NamedObject dateObj = 
(NamedObject)grid.getTimes().get(0);
                                                 String dateStr = 
dateObj.getDescription();
                                                 // dateStr should be formatted 
to 20 chars
                                                 // time format in our netCDF 
files 
                                                 // days_since_20 10101.06
                                                 // format in tiff is 
YYYY:MM:DD HH:MM:SS
                                                 dateStr = dateStr.replace(' ', 
'0');
                                                 DateFormat formatter = new 
SimpleDateFormat("'days_since_'yyyyMMdd'.'HH");
                                                 Date date = 
formatter.parse(dateStr);
                                                 DateFormat tiffFormatter = new 
SimpleDateFormat("yyyy:MM:dd HH:mm:ss");
                                                 // null terminated string
                                                 String[] tiffDateStr = 
{tiffFormatter.format(date)};
                                                 dateField = new 
TIFFField(BaselineTIFFTagSet.TAG_DATE_TIME, TIFFField.TIFF_ASCII, 1, 
tiffDateStr);
                                                 }*/
                                                if (dateField != null) {
                                                        TIFFField[] fields = { 
descField, geoKeyDir,
                                                                        
tiePtField, pixelSclField, dateField };
                                                        
params.setExtraFields(fields);
                                                } else {
                                                        TIFFField[] fields = { 
descField, geoKeyDir,
                                                                        
tiePtField, pixelSclField };
                                                        
params.setExtraFields(fields);
                                                }

                                        }

                                        
params.setCompression(TIFFEncodeParam.COMPRESSION_NONE);
                                        ImageEncoder encoder = 
ImageCodec.createImageEncoder(
                                                        "TIFF", out, params);

                                        if (encoder == null) {
                                                
System.out.println("imageEncoder is null");
                                                System.exit(0);
                                        }
                                        encoder.encode(result);
                                        out.close();
                                }
                        }  // end of while
                        if (!found)
                        {
                                System.out.println("failed to find named 
variable in dataset");
                        }
                } catch (Exception e) {
                        e.printStackTrace();
                }

        }
        
        public static void printUsage()
        {
                System.out.println("Usage: java TestGeotiffWriter netCDFFile 
variableName");
        }
        
        public static void main(String[] args)
        {
                if (args.length < 2)
                {
                        TestGeotiffWriter.printUsage();
                }
                else
                {
                        TestGeotiffWriter writer = new TestGeotiffWriter();
                        writer.convertToGeotiff(args[0], args[1]);
                        System.out.println("Finished conversion");
                }
        }

}


-----Original Message-----
From: John Caron [mailto:address@hidden]
Sent: Tuesday, April 19, 2005 3:13 PM
To: Norman Barker
Cc: Yuan Ho; address@hidden;
address@hidden
Subject: Re: [Fwd: multi band tiff from netCDF]


Hi Norman:

Thanks for the code. We looked at using JAI to write TIFF, but failed to 
get it working ("JAI documentation is a bit sparse"). Is the code the 
complete way to do that? Are there any problems with it that you know 
of? It doesnt appear that youve added the geotiff tags, is that correct? 
Any advice you can give us on using or not using JAI?

thanks again!




Norman Barker wrote:

>Thanks for the reply, I will top post (apologies) since the file is long.
>
>I wrote some test code before I integrated into my app, and using netcdf-java 
>and
>JAI the following creates a multiband tif file. (there are some redundant 
>variables in
>there, apologies again).  JAI documentation is a bit sparse, so hope this 
>helps.
>
>Norman
>
>       public static void main(String[] args) {
>               String testNC = "/temp/test.nc";
>               String outGTif = "/temp/out.tif";
>               String gridName = "ta";
>               try {
>                       NetcdfDataset ncDataset = 
> NetcdfDataset.openDataset(testNC);
>                       
>                       GridDataset gSet = new GridDataset(ncDataset);
>                       GeoGrid grid = gSet.findGridByName(gridName);
>            GridCoordSys gcs = grid.getCoordinateSystem();
>                       
>                   //latlon coord does not need to be scaled
>                   double scaler = (gcs.isLatLon()) ? 1.0 : 1000.0;
>
>                   CoordinateAxis1D xaxis = (CoordinateAxis1D) 
> gcs.getXHorizAxis();
>                   CoordinateAxis1D yaxis = (CoordinateAxis1D) 
> gcs.getYHorizAxis();
>                   double xStart = xaxis.getCoordValue(0) * scaler;
>                   double yStart = yaxis.getCoordValue(0) * scaler;
>                       double xEnd = 
> xaxis.getCoordValue(grid.getXDimension().getLength()-1) * scaler;
>                       double yEnd = 
> yaxis.getCoordValue(grid.getYDimension().getLength()-1) * scaler;
>                       double xInc = (xEnd - xStart)/xaxis.getElementSize();
>                       double yInc = (yEnd - yStart)/yaxis.getElementSize();
>                       
>                       
>                       
>                       List list = ncDataset.getCoordinateAxes();
>                       Iterator itr = list.iterator();
>
>                       while (itr.hasNext()) {
>                               CoordinateAxis1D dim = (CoordinateAxis1D) 
> itr.next();
>                               String dimName = dim.getName();
>
>                               if (dimName.equals(IngestFile.DIM_NAME)) {
>                                       int nBands = (int) dim.getSize();
>                                       if (dim.isUnlimited()) {
>                                               // we can't do anything, log
>                                               return;
>                                       }
>                                       
>                                       int width = 
> grid.getXDimension().getLength();
>                                       int height = 
> grid.getYDimension().getLength();
>               
>                                       ParameterBlockJAI pbjai = new 
> ParameterBlockJAI("bandmerge");
>                                       
>                                       for (int i = 0; i < nBands; i++)
>                                       {
>                                               Array arr = 
> grid.readDataSlice(-1, i, -1, -1);
>                                               float[] data = (float[]) 
> arr.getStorage();
>                                               DataBufferFloat dbuffer = new 
> DataBufferFloat(data, width*height);
>                                               
>                                               // create a float data sample 
> model
>                                               SampleModel sampleModel = 
> RasterFactory.createBandedSampleModel(DataBuffer.TYPE_FLOAT, width, height, 
> 1);
>                                               
>                                               // create a compatible colour 
> model
>                                               ColorModel colorModel = 
> PlanarImage.createColorModel(sampleModel);
>                                               
>                                               // create a writable raster
>                                               Raster raster = 
> RasterFactory.createWritableRaster(sampleModel, dbuffer, null);
>                                               
>                                               // create a tiled image using 
> the float sample model
>                                               TiledImage tiledImage = new 
> TiledImage(0, 0, width, height, 0, 0, sampleModel, colorModel);
>                                               
>                                               // set the data of the tiled 
> image to be the raster
>                                               tiledImage.setData(raster);
>                                               
>                                               pbjai.setSource(tiledImage, i);
>                                       }
>                                       // save the image on a file
>                                       PlanarImage result = 
> JAI.create("bandmerge", pbjai, null);
>                                       JAI.create("filestore", result, 
> outGTif, "TIFF");
>                                       break;
>                               }
>                       }
>               } catch (Exception e) {
>                       e.printStackTrace();
>               }
>
>       }
>
>-----Original Message-----
>From: Yuan Ho [mailto:address@hidden]
>Sent: Monday, April 18, 2005 5:27 PM
>To: Norman Barker
>Cc: address@hidden; address@hidden
>Subject: Re: [Fwd: multi band tiff from netCDF]
>
>
>John Caron wrote:
>
>  
>
>>can you answer this?
>>cc to address@hidden and 
>>address@hidden
>>
>>-------- Original Message --------
>>Subject:     multi band tiff from netCDF
>>Date:     Mon, 18 Apr 2005 16:01:21 +0100
>>From:     Norman Barker <address@hidden>
>>Organization:     UCAR/Unidata
>>To:     <address@hidden>
>>
>>
>>
>>Hi,
>>
>>I am using the java classes  to convert a netcdf file for one time 
>>instance, with 27 levels in Z (atmospheric pressure levels).
>>
>>I have created a loop with variable i, and attempting to write a band 
>>in the tif file
>>for each iteration. What seems to happen is the first band gets 
>>written ok, but then after that
>>the page number gets messed up?
>>
>>Array data = grid.readDataSlice(0, i, -1, -1);
>>writer.writeGrid(grid, data, false, xStart, yStart, xInc, yInc, i);
>>
>>Has anyone got any pseudo code for writing multiband tif files from 
>>netCDF?
>>
>>Many thanks,
>>
>>Norman Barker
>>    
>>
>
>Norman,
>               The code was not designed for multiband tiff files. We 
>originally believed that multipage was better to store multi variables 
>and/or multi vertical levels.
>We just realized that multiband was probably better embraced because of 
>the existing tiff readers and the concept of tiff community. This will 
>be my next task to change the
>software to use multiband to store tiff files, it will be released 
>sometime in the Summer.
>
>Yuan
>
>  
>