Delbert Franz wrote:
> I have been working to display an image of a USGS 7.5 minute quad sheet.
> These are provided at various locations about the Web. Since the
> range of colors on these maps is limited, the *.tif files appear to
> use an indexed color map wherein each pixel has a value 0 to 255 and
> the color is found from a table with 256 entries having triplets of
> the RGB in the range of 0-255. I have not been able to sort out
> how to get the gdal package to give me the color map from within python, so I
> dumped it
> from the image file using gdalinfo and then cut and pasted to get
> the following script:
>
> ---------------------------------------------------------------------------------
> import numpy as np
> import matplotlib.pyplot as plt
> import matplotlib.cm as cm
> from matplotlib.colors import ListedColormap
from matplotlib.colors import NoNorm
>
> import osgeo.gdal as gdal
> from osgeo.gdalconst import *
>
> gd = gdal.Open('o37122d1.tif')
>
> #Setup to compute the colormap from the RGB triplets from the geotif file
> div = np.zeros( (256,3), np.float32)
>
> div = 255.0
> ctab = np.array([
> [ 255,255,255],
> [ 0,0,0],
> [ 255,255,255],
> [ 91,159,230],
> [ 230,45,30],
> [ 162,96,71],
> [ 210,255,177],
> [ 197,101,197],
> [ 255,240,0],
> [ 202,225,245],
>
> ...snip.... (deleted many lines:)
> [ 250,202,250],
> [ 230,230,230],
> [ 222,167,146],
> [ 255,255,255],
> [ 255,255,255] ], dtype=np.uint8)
>
> #Compute colors in range 0.0 to 1.0.
> fctab= ctab/div
>
> usgscm = ListedColormap(fctab, name='usgs',N=None)
>
>
> doq = gd.ReadAsArray()
> # doq proves to be a uint8 array: 8802 rows and 7066 columns
>
> # Cut out a subset from the main array--to large to display and
> # slow as 'molasses in January' to process:)
> suba = doq[0:1101 ,0:1801 ]
>
> fig = plt.figure()
> ax = fig.add_subplot(111)
> ax.imshow(suba, cmap=usgscm, origin='upper')
Instead, try:
ax.imshow(suba, cmap=usgscm, norm=NoNorm(), origin='upper')
>
> plt.show()
> ---------------------------------------------------------------------------------
> This script does give me an image--but in badly wrong colors:( The script
> does
> properly display gray-scaled Digital Ortho-quadrangles using cm.gray as the
> color
> map in imshow. Consequently something is not quite correct with respect to
> the
> definition or the use of the color map. It appears that each map, and there
> are about 56,000 of them available on one site, could have its own color map.
> Thus my application must be able to compute a color map unique to each of the
> topographic maps.
>
> Questions:
>
> 1. What am I missing to get imshow to pick out the correct colors from the
> color map?
The default norm will scale your inputs; specifying the norm as a NoNorm
instance will pass the integers through directly, so they will be used
as indices into the colormap.
>
> 2. Should I be using the name, usgs, given in the ListedColormap instance
> someplace?
> Not clear to me what role that name plays.
>
None, really. I can imagine ways in which it could be useful, but
unless you know you need it, consider it optional.
> 3. Is there a way to use ctab directly in the ListedColormap instance? Class
> Colormap
> has a bytes argument, False by default, but I am not yet sure if it has any
> bearing
> on my problem.
No, sorry, but the bytes argument is only in the __call__ method, where
it is used to improve efficiency within mpl. There is no facility for
using ints in the lookup table, and no recognition of 0-255 ints within
mpl colorspecs, hence no way to feed your ctab in directly. Even if
there were, though, I don't think it would make much difference in
plotting time.
Eric
>
> I am using Matplotlib 0.98 with the latest numpy, and gdal. I am not using
> Basemap, after some study, because it contains far more "horsepower" than my
> simple topographic-map viewer application requires. Also, this is my first
> foray into "image processing". I am finding the learning curve a bit steep,
> as usual, with multiple names for the same thing popping up in descriptions
> from various sources--business as usual in the computer business:)
>
> Thanks from a happy matplotlib user in California!
>
> Delbert
>
>
> -------------------------------------------------------------------------
> Check out the new SourceForge.net Marketplace.
> It's the best place to buy or sell services for
> just about anything Open Source.
> http://sourceforge.net/services/buy/index.php
> _______________________________________________
> Matplotlib-users mailing list
> Matplotlib-users@lists.sourceforge.net
> https://lists.sourceforge.net/lists/listinfo/matplotlib-users
-------------------------------------------------------------------------
Check out the new SourceForge.net Marketplace.
It's the best place to buy or sell services for
just about anything Open Source.
http://sourceforge.net/services/buy/index.php
_______________________________________________
Matplotlib-users mailing list
Matplotlib-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/matplotlib-users
