AllSkyVisualization: makeWebGL.py

File makeWebGL.py, 2.0 KB (added by Jayce Dowell, 7 years ago)

Script to generate a WebGL base image from a HEALpix map

Line 
1#!/usr/bin/env python
2
3import os
4import sys
5import numpy
6import healpy
7
8def imsave(fname, arr, vmin=None, vmax=None, cmap=None, format=None, origin=None):
9    from matplotlib.backends.backend_agg import FigureCanvasAgg as FigureCanvas
10    from matplotlib.figure import Figure
11
12    fig = Figure(figsize=arr.shape[::-1], dpi=1, frameon=False)
13    canvas = FigureCanvas(fig)
14    fig.figimage(arr, cmap=cmap, vmin=vmin, vmax=vmax, origin=origin)
15    fig.savefig(fname, dpi=1, format=format)
16
17from matplotlib.mlab import griddata
18import matplotlib.pyplot as plt
19import matplotlib.cm as cm
20from matplotlib.ticker import NullFormatter
21
22import healpy.rotator as R
23from healpy import pixelfunc
24from healpy.pixelfunc import UNSEEN
25
26
27def main(args):
28        filename = args[0]
29        map = numpy.loadtxt(filename)
30        wgt = map*0.0 + 1.0
31       
32        filename = filename.replace('map', 'wgt')
33        if os.path.exists(filename):
34                wgt = numpy.loadtxt(filename)
35                map /= wgt
36
37        map = healpy.pixelfunc.ud_grade(map, 128)
38        wgt= healpy.pixelfunc.ud_grade(wgt, 128)
39        mask = numpy.where( wgt == 0, 1, 0 )
40       
41        map = healpy.ma(map)
42        wgt = healpy.ma(wgt)
43        map.mask = mask
44        wgt.mask = mask
45       
46        good = numpy.ma.array(map.filled(), mask=mask)
47        map = healpy.smoothing(map.filled(map.mean()), fwhm=1.0, degree=True)
48       
49        map = healpy.ma(map)
50        map.mask = mask
51
52        # Basename
53        basename = os.path.split(filename)[1]
54        basename = os.path.splitext(basename)[0]
55        j1, j2, j3, freq1, freq2 = basename.split('-', 4)
56        print freq1, freq2
57
58        # Compute a reasonable color map to display the image with
59        from scipy.stats import scoreatpercentile as percentile
60        vmin = percentile(good.compressed(), 10)
61        vmax = percentile(good.compressed(), 99.5)
62
63        # Project and save
64        p = healpy.projector.CartesianProj(xsize=1024)
65        nside = healpy.pixelfunc.npix2nside(healpy.pixelfunc.get_map_size(map))
66        f = lambda x,y,z: healpy.pixelfunc.vec2pix(nside,x,y,z)
67        imsave("lwa%s.png" % freq1, p.projmap(map, f), vmin=vmin, vmax=vmax, origin='lower')
68
69       
70if __name__ == "__main__":
71        main(sys.argv[1:])
72