Using indexed color maps for images-not getting correct colors

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

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')

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?

2. Should I be using the name, usgs, given in the ListedColormap instance someplace?
   Not clear to me what role that name plays.

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.

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

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

Eric,

That did the trick! The images are displaying properly now. Thanks
for the quick help.

                             Delbert

···

On Monday 23 June 2008 23:00, Eric Firing wrote:

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)
>
> 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

>