- #!/usr/bin/env python
- """Generate an image from ASCII data.
- Step 1: make this file executable::
- chmod +x plot_image.py
- Step 2: pipe data into python script::
- ./gen_data | ./plot_image.py -s nx,ny -c nc -t 'image title'
- Data should be one ASCII float per line, in the following order::
- z[x1,y1]
- z[x2,y1]
- ...
- z[xN,y1]
- z[x1,y2]
- ...
- z[xN,y2]
- ...
- z[xN,yN]
- where `x` increases from `x1` to `xN` and `y` increases from `y1`
- through `yN`.
- You can use the `--xyz` option to read an alternative data format,
- which is::
- x1 y1 z[x1,y1]
- x2 y2 z[x1,y1]
- ...
- If you use this method, the ordering of lines in the data file is
- irrelevant, and `Nx` and `Ny` are extracted from the data file itself.
- However, you still need to use a rectangular grid (i.e. for every
- `xi`, you need to have entries for every `yi`).
- For usage info, run::
- ./plot_image.py --help
- When run with interactive output (i.e. no `-o ...` option), the
- interactive figure is displayed using `pylab.show`, which means that
- you'll have to kill `plot_image.py` using ^C or similar [1]_.
- For other ideas, see the Matplotlib website [2]_.
- .. [1] http://matplotlib.sourceforge.net/faq/howto_faq.html#use-show
- .. [2] http://matplotlib.sourceforge.net/
- """
- import optparse
- import sys
- import matplotlib
- import matplotlib.image
- import numpy
- # Depending on your Matplotlib configuration, you may need to adjust
- # your backend. Do this before importing pylab or matplotlib.backends.
- #matplotlib.use('Agg') # select backend that doesn't require X Windows
- #matplotlib.use('GTKAgg') # select backend that supports pylab.show()
- import pylab
- _DOC = __doc__
- def read_data_1d(stream, nx, ny):
- """Read in data, one entry per line.
- >>> from StringIO import StringIO
- >>> s = StringIO('\\n'.join(map(str, range(10)))+'\\n')
- >>> X,Y,Z = read_data_1d(s, 5, 2)
- >>> X
- array([[0, 1, 2, 3, 4, 5],
- [0, 1, 2, 3, 4, 5],
- [0, 1, 2, 3, 4, 5]])
- >>> Y
- array([[0, 0, 0, 0, 0, 0],
- [1, 1, 1, 1, 1, 1],
- [2, 2, 2, 2, 2, 2]])
- >>> Z
- array([[ 0., 1., 2., 3., 4.],
- [ 5., 6., 7., 8., 9.]])
- """
- X,Y = pylab.meshgrid(range(nx+1), range(ny+1))
- Z = numpy.loadtxt(stream)
- assert Z.size == nx*ny, 'Z.size = %d != %d = %dx%d' % (
- Z.size, nx*ny, nx, ny)
- Z = Z.reshape([x-1 for x in X.shape])
- return (X,Y,Z)
- def read_data_3d(stream):
- """Read in data, one `(x, y, z)` tuple per line.
- >>> from StringIO import StringIO
- >>> lines = []
- >>> for x in range(5):
- ... for y in range(2):
- ... lines.append('\t'.join(map(str, [x, y, x+y*5])))
- >>> s = StringIO('\\n'.join(lines)+'\\n')
- >>> X,Y,Z = read_data_3d(s)
- >>> X
- array([[ 0., 1., 2., 3., 4., 5.],
- [ 0., 1., 2., 3., 4., 5.],
- [ 0., 1., 2., 3., 4., 5.]])
- >>> Y
- array([[ 0., 0., 0., 0., 0., 0.],
- [ 1., 1., 1., 1., 1., 1.],
- [ 2., 2., 2., 2., 2., 2.]])
- >>> Z
- array([[ 0., 1., 2., 3., 4.],
- [ 5., 6., 7., 8., 9.]])
- """
- XYZ = numpy.loadtxt(stream)
- assert len(XYZ.shape) == 2 and XYZ.shape[1] == 3, XYZ.shape
- Xs = numpy.array(sorted(set(XYZ[:,0])))
- Ys = numpy.array(sorted(set(XYZ[:,1])))
- Z = numpy.ndarray((len(Ys), len(Xs)), dtype=float)
- xyz = {} # dict of z values keyed by (x,y) tuples
- for i in range(XYZ.shape[0]):
- xyz[(XYZ[i,0], XYZ[i,1])] = XYZ[i,2]
- for i,x in enumerate(Xs):
- for j,y in enumerate(Ys):
- Z[j,i] = xyz[x,y]
- # add dummy row/column for pcolor
- dx = Xs[-1] - Xs[-2]
- dy = Ys[-1] - Ys[-2]
- Xs = numpy.append(Xs, Xs[-1] + dx)
- Ys = numpy.append(Ys, Ys[-1] + dy)
- X,Y = pylab.meshgrid(Xs, Ys)
- return (X,Y,Z)
- def plot(X, Y, Z, full=False, title=None, contours=None, interpolation=None,
- cmap=None):
- """Plot Z over the mesh X, Y.
- >>> X, Y = pylab.meshgrid(range(6), range(2))
- >>> Z = X[:-1,:-1]**2 + Y[:-1,:-1]
- >>> plot(X, Y, Z) # doctest: +ELLIPSIS
- <matplotlib.figure.Figure object at 0x...>
- """
- X_min = X[0,0]
- X_max = X[-1,-1]
- Y_min = Y[0,0]
- Y_max = Y[-1,-1]
- fig = pylab.figure()
- if full:
- axes = fig.add_axes([0, 0, 1, 1])
- else:
- axes = fig.add_subplot(1, 1, 1)
- if title:
- axes.set_title(title)
- axes.set_axis_off()
- if contours:
- cset = axes.contour(X[:-1,:-1], Y[:-1,:-1], Z, contours, cmap=cmap)
- # [:-1,:-1] to strip dummy last row & column from X&Y.
- axes.clabel(cset, inline=1, fmt='%1.1f', fontsize=10)
- else:
- # pcolor() is much slower than imshow.
- #plot = axes.pcolor(X, Y, Z, cmap=cmap, edgecolors='none')
- #axes.autoscale_view(tight=True)
- plot = axes.imshow(Z, aspect='auto', interpolation=interpolation,
- origin='lower', cmap=cmap,
- extent=(X_min, X_max, Y_min, Y_max))
- if not full:
- fig.colorbar(plot)
- return fig
- def get_possible_interpolations():
- try: # Matplotlib v1.0.1
- return sorted(matplotlib.image.AxesImage._interpd.keys())
- except AttributeError:
- try: # Matplotlib v0.91.2
- return sorted(matplotlib.image.AxesImage(None)._interpd.keys())
- except AttributeError:
- # give up ;)
- pass
- return ['nearest']
- def test():
- import doctest
- results = doctest.testmod()
- return results.failed
- def main(argv=None):
- """Read in data and plot it.
- >>> from tempfile import NamedTemporaryFile
- >>> i = NamedTemporaryFile(prefix='tmp-input', suffix='.dat')
- >>> i.write('\\n'.join([str(x) for x in range(10)])+'\\n')
- >>> i.flush()
- >>> o = NamedTemporaryFile(prefix='tmp-output', suffix='.png')
- >>> main(['-i', i.name, '-s', '5,2', '-o', o.name, '-m', 'binary'])
- Plot_image
- Title: Some like it hot
- Image size: 5 2
- False color
- X range: 0 4 (6 steps)
- Y range: 0 1 (3 steps)
- Z range: 0.0 9.0
- >>> img = o.read()
- >>> img.startswith('\\x89PNG')
- True
- >>> i.close()
- >>> o.close()
- """
- if argv == None:
- argv = sys.argv[1:]
- usage = '%prog [options]'
- epilog = _DOC
- p = optparse.OptionParser(usage=usage, epilog=epilog)
- p.format_epilog = lambda formatter: epilog+'\n'
- p.add_option('-i', '--input', dest='input', metavar='FILE',
- help='If set, read data from FILE rather than stdin.')
- p.add_option('-o', '--output', dest='output', metavar='FILE',
- help=('If set, save the figure to FILE rather than '
- 'displaying it immediately'))
- p.add_option('-s', '--size', dest='size', default='%d,%d' % (16, 16),
- help='Data size (columns,rows; default: %default)')
- p.add_option('-3', '--xyz', dest='xyz', default=False, action='store_true',
- help=('If set, read (x,y,z) tuples from the input data rather'
- 'then reading `z` and calculating `x` and `y` from '
- '`--size`.'))
- p.add_option('-c', '--contours', dest='contours', type='int',
- help=('Number of contour lines (if not set, draw false color '
- 'instead of contour lines; default: %default)'))
- p.add_option('-f', '--full-figure', dest='full', action='store_true',
- help=('Set axes to fill the figure (i.e. no title or color '
- 'bar'))
- p.add_option('-t', '--title', dest='title', default='Some like it hot',
- help='Title (%default)')
- p.add_option('--test', dest='test', action='store_true',
- help='Run internal tests and exit.')
- interpolations = get_possible_interpolations()
- p.add_option('--interpolation', dest='interpolation', default='nearest',
- help=('Interpolation scheme (for false color images) from %s '
- '(%%default)') % ', '.join(interpolations))
- maps=[m for m in pylab.cm.datad if not m.endswith("_r")]
- maps.sort()
- p.add_option('-m', '--color-map', dest='cmap', default='jet',
- help='Select color map from %s (%%default)' % ', '.join(maps))
- options,args = p.parse_args(argv)
- if options.test:
- sys.exit(test())
- nx,ny = [int(x) for x in options.size.split(',')]
- try:
- cmap = getattr(pylab.cm, options.cmap)
- except AttributeError:
- raise Exception('no color map named %s in %s'
- % (options.cmap, ', '.join(maps)))
- print 'Plot_image'
- print 'Title: ', options.title
- if not options.xyz:
- print 'Image size: ', nx, ny
- if options.contours:
- print '# countour lines: ', options.contours
- else:
- print 'False color'
- if options.input:
- fin = open(options.input, 'r')
- else:
- fin = sys.stdin
- if options.xyz:
- X,Y,Z = read_data_3d(fin)
- else:
- X,Y,Z = read_data_1d(fin, nx, ny)
- if options.input:
- fin.close()
- Z_min = numpy.min(Z.flat)
- Z_max = numpy.max(Z.flat)
- print 'X range: {} {} ({} steps)'.format(
- X[0,0], X[0,-2], X.shape[1])
- print 'Y range: {} {} ({} steps)'.format(
- Y[0,0], Y[-2,0], Y.shape[0])
- print 'Z range: ', Z_min, Z_max
- fig = plot(X, Y, Z, full=options.full, title=options.title,
- contours=options.contours, interpolation=options.interpolation,
- cmap=cmap)
- if options.output:
- fig.savefig(options.output)
- else:
- pylab.show()
- if __name__ == '__main__':
- main()