<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">
<HTML><HEAD>
<META content="text/html; charset=utf-8" http-equiv=Content-Type>
<META name=GENERATOR content="MSHTML 8.00.6001.18928"></HEAD>
<BODY>
<DIV dir=ltr align=left><SPAN class=441562615-22092010><FONT size=2
face=Arial>Chaitanya,</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=441562615-22092010><FONT size=2
face=Arial>Thanks for the help, but it didn't change anything. If it helps, I'm
using gdal 1.7.2 with python 3.1 and numpy 1.5.</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=441562615-22092010><FONT size=2
face=Arial>Cheers,</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=441562615-22092010><FONT size=2
face=Arial>Andy</FONT></SPAN></DIV><BR>
<DIV dir=ltr lang=en-us class=OutlookMessageHeader align=left>
<HR tabIndex=-1>
<FONT size=2 face=Tahoma><B>From:</B> Chaitanya kumar CH
[mailto:chaitanya.ch@gmail.com] <BR><B>Sent:</B> 22 September 2010
16:11<BR><B>To:</B> Hartley, Andrew<BR><B>Cc:</B>
gdal-dev@lists.osgeo.org<BR><B>Subject:</B> Re: [gdal-dev] RE: Compression using
the create method in python and aggregation methods<BR></FONT><BR></DIV>
<DIV></DIV>Andy,<BR><BR>Try this instead.<BR><FONT
style="FONT-FAMILY: courier new,monospace" size=2>out =
outDrv.Create("outfile.tif", gscl.RasterXSize, gscl.RasterYSize, 1,
gdalconst.GDT_UInt16, </FONT><B
style="FONT-FAMILY: courier new,monospace"><FONT size=2>options =
</FONT></B><FONT style="FONT-FAMILY: courier new,monospace" size=2>[
'COMPRESS=LZW' ] )</FONT> <BR><BR>However, I am not sure this makes any
difference because gdalinfo did tell that there is LZW compression.<BR><BR>
<DIV class=gmail_quote>On Wed, Sep 22, 2010 at 4:41 PM, Hartley, Andrew <SPAN
dir=ltr><<A
href="mailto:andrew.hartley@metoffice.gov.uk">andrew.hartley@metoffice.gov.uk</A>></SPAN>
wrote:<BR>
<BLOCKQUOTE
style="BORDER-LEFT: rgb(204,204,204) 1px solid; MARGIN: 0pt 0pt 0pt 0.8ex; PADDING-LEFT: 1ex"
class=gmail_quote>
<DIV>
<P><FONT size=2 face=Arial>Hi all,</FONT> </P>
<P><FONT size=2 face=Arial>I'm trying to create a Gtiff with LZW compression
using python, with the code below, which I wrote with help from the tutorial
at </FONT><A href="http://www.gdal.org/gdal_tutorial.html"
target=_blank><U><FONT color=#0000ff size=2
face=Arial>http://www.gdal.org/gdal_tutorial.html</FONT></U></A><FONT size=2
face=Arial>. Gdalinfo tells me that the resulting tif ("outfile.tif") has
compression (Image Structure Metadata: COMPRESSION=LZW), but all my outfiles
have the same file size and are quite large, so it seems they actually aren't
compressed. In fact, when I tried:</FONT></P>
<P><FONT size=2 face=Arial>gdal_translate -co 'COMPRESS=LZW' outfile.tif
newoutfile.tif</FONT> </P>
<P><FONT size=2 face=Arial>the newoutfile.tif is considerably smaller than the
original file. So, this leads me to think that there's a problem with my
create() statement below (see the line in bold below). Could somebody please
tell me what I have missed?</FONT></P>
<P><FONT size=2 face=Arial>Since I'm here, I think I will also pick your
brains about cell aggregation methods. You'll see from the code below that I
have writen a loop to aggregate only cells with data. I spent a bit of time
considering a few options (for example, the excellent pages by Dr Gomez-Dans -
</FONT><A
href="http://sites.google.com/site/spatialpython/aggregating-data-to-grid-cells"
target=_blank><U><FONT color=#0000ff size=2
face=Arial>http://sites.google.com/site/spatialpython/aggregating-data-to-grid-cells</FONT></U></A><FONT
size=2 face=Arial>). My code works reasonably well, but since I have lots of
processing to do and it is not as fast as I would like, I was wondering if
anybody could suggest a more efficient solution?</FONT></P>
<P><FONT size=2 face=Arial>Thanks very much in advance for any help you may be
able to offer me!</FONT> </P>
<P><FONT size=2 face=Arial>Kind regards,</FONT> <BR><FONT size=2
face=Arial>Andy</FONT> </P>
<P><FONT size=2 face=Arial>s = (640,640)</FONT> <BR><FONT size=2 face=Arial>dt
= numpy.dtype('uint16') </FONT><BR><FONT size=2
face=Arial># reftile is approx 1km resolution raster, with a unique ID for
each cell for a 5 degree square window</FONT> <BR><FONT size=2 face=Arial>gscl
= gdal.Open (reftile)</FONT> <BR><FONT size=2 face=Arial>tilescl =
g.GetRasterBand(1).ReadAsArray().astype(numpy.uint16) </FONT><BR><FONT size=2
face=Arial># reftile90 is a 90m resample (using nearest neighbour) of
reftile</FONT> <BR><FONT size=2 face=Arial>g90 = gdal.Open (reftile90)</FONT>
<BR><FONT size=2 face=Arial>tile90 =
g.GetRasterBand(1).ReadAsArray().astype(numpy.uint16)
</FONT><BR>
<BR><FONT size=2 face=Arial>z =
numpy.zeros(s, dtype=dt)</FONT> <BR><FONT size=2 face=Arial>U =
unique(tile90[numpy.greater(rec90, 0)])</FONT> <BR><FONT size=2
face=Arial>lenU = len(U)</FONT> <BR><FONT size=2 face=Arial># for each 90m
cell with data, aggregate and write to low resolution output grid</FONT>
<BR><FONT size=2 face=Arial>for u in range(lenU):</FONT>
<BR> <FONT size=2 face=Arial>result
= numpy.sum(rec90[numpy.equal(tile90, U[u])])</FONT>
<BR> <FONT size=2
face=Arial>z[numpy.equal(tilescl, U[u])] = result</FONT> <BR><FONT size=2
face=Arial># Write out the grid</FONT> <BR><FONT size=2 face=Arial>outDrv =
gdal.GetDriverByName('GTiff')</FONT> <BR><B><FONT size=2 face=Arial>out =
outDrv.Create("outfile.tif", gscl.RasterXSize, gscl.RasterYSize, 1,
gdalconst.GDT_UInt16, [ 'COMPRESS=LZW' ] )</FONT></B> <BR><FONT size=2
face=Arial>out.SetProjection(gscl.GetProjection())</FONT> <BR><FONT size=2
face=Arial>out.SetGeoTransform(gscl.GetGeoTransform())</FONT> <BR><FONT size=2
face=Arial>out.GetRasterBand(1).WriteArray(z)</FONT> <BR><FONT size=2
face=Arial>gscl = None</FONT> <BR><FONT size=2 face=Arial>G90 = None</FONT>
<BR><FONT size=2 face=Arial>out = None</FONT> </P><BR><FONT color=#888888><BR>
<P><FONT size=2 face=Arial>--</FONT> <BR><FONT size=2 face=Arial>Andrew
Hartley Climate Impacts Risk Analyst</FONT> <BR><B><FONT size=2
face=Arial>Met Office Hadley Centre</FONT></B><FONT size=2 face=Arial>
FitzRoy Road Exeter Devon EX1 3PB United
Kingdom</FONT> <BR><FONT size=2 face=Arial>Tel: +44 (0)1392 885720 Fax:
+44 (0)1392 885681</FONT> <BR><FONT size=2 face=Arial>Email: <A
href="mailto:andrew.hartley@metoffice.gov.uk"
target=_blank>andrew.hartley@metoffice.gov.uk</A> Website:
</FONT><A><U><FONT color=#0000ff size=2
face=Arial>www.metoffice.gov.uk</FONT></U></A> </P>
<P><FONT size=2 face=Arial>See our guide to climate change at </FONT><A
href="http://www.metoffice.gov.uk/climatechange/guide/" target=_blank><U><FONT
color=#0000ff size=2
face=Arial>http://www.metoffice.gov.uk/climatechange/guide/</FONT></U></A>
</P></FONT></DIV><BR>_______________________________________________<BR>gdal-dev
mailing list<BR><A
href="mailto:gdal-dev@lists.osgeo.org">gdal-dev@lists.osgeo.org</A><BR><A
href="http://lists.osgeo.org/mailman/listinfo/gdal-dev"
target=_blank>http://lists.osgeo.org/mailman/listinfo/gdal-dev</A><BR></BLOCKQUOTE></DIV><BR><BR
clear=all><BR>-- <BR>Best regards,<BR>Chaitanya kumar CH.<BR>/tʃaɪθənjə/ /kʊmɑr/
<BR>+91-9494447584<BR>17.2416N 80.1426E<BR></BODY></HTML>