working_with_wham.ipynb
from matplotlib import pyplot as plt
import numpy as np
from astropy.io import fits
import healpy as hp
from healpy.newvisufunc import projview, newprojplot
from astropy.wcs import WCS
from astropy.coordinates import SkyCoord
import astropy.units as u
import h5py
from scipy.ndimage import map_coordinates
fn = '/home/justin/Research/WHAM/wham-ss-DR1-v161116-170912-grid.fits'
whamdat = fits.getdata(fn)
whamhdr = fits.getheader(fn)
whamhdr
WARNING: VerifyWarning: Verification reported errors: [astropy.io.fits.verify]
WARNING: VerifyWarning: Card 'RESTWAV' is not FITS standard (invalid value string: '6.56300e-07 / Rest wavelength of H-alpha [m]'). Fixed 'RESTWAV' card to meet the FITS standard. [astropy.io.fits.verify]
WARNING: VerifyWarning: Note: astropy.io.fits uses zero-based indexing.
[astropy.io.fits.verify]
SIMPLE = T / Written by IDL: Tue Sep 12 13:28:51 2017
BITPIX = -32 / Real*4 (floating point)
NAXIS = 3 /
NAXIS1 = 1440 /
NAXIS2 = 720 /
NAXIS3 = 162 /
EXTEND = T / Extensions may be present
DATE = '2017-09-12' /
TAG = 'DR1 ' / General release tag
VERSION = '161116 ' / Generated from wham-ss-161116.sav
OBJECT = 'WHAM Sky Survey - Release DR1 (161116)' /
RESTWAV = 6.56300E-07 / Rest wavelength of H-alpha [m]
BSCALE = 1.00000 /
BZERO = 0.00000 /
BUNIT = 'mR (m s^-1)^-1' / mili-Rayleighs per m/s, equiv. R per km/s
CTYPE1 = 'GLON-CAR' / Coordinate Type
CTYPE2 = 'GLAT-CAR' / Coordinate Type
EQUINOX = 2000.00000000 /Equinox of Ref. Coord.
CD1_1 = -0.250000 / Degrees / Pixel
CD2_1 = 0.00000 / Degrees / Pixel
CD1_2 = 0.00000 / Degrees / Pixel
CD2_2 = 0.250000 / Degrees / Pixel
CRPIX1 = 720.500 / Reference Pixel in X
CRPIX2 = 360.500 / Reference Pixel in Y
CRVAL1 = 0.00000000000 / Galactic longitude of reference pixel
CRVAL2 = 0.00000000000 / Galactic latitude of reference pixel
PV1_0 = 0.00000000000 / Projection parameters
PV1_1 = 0.00000000000 / Projection parameters
PV1_2 = 0.00000000000 / Projection parameters
PV1_3 = 0.00000000000 / Projection parameters
PV1_4 = 90.0000000000 / Projection parameters
CTYPE3 = 'VELO-LSR' /
SPECSYS = 'LSRK ' /
CRVAL3 = -169116. /
CRPIX3 = 1.00000 /
CDELT3 = 2027.58 /
CUNIT1 = 'deg ' /
CUNIT2 = 'deg ' /
CUNIT3 = 'm/s ' /
CROTA3 = 0.00000 /
COMMENT WHAM Sky Survey - Release DR1 (161116)
COMMENT
COMMENT REFERENCES:
COMMENT
COMMENT Northern Sky Survey - Haffner, Reynolds, Tufte, Madsen, Jaehnig, &
COMMENT Percival. 2003, ApJS, 149, 405, doi:10.1086/378850.
COMMENT
COMMENT (Full) Sky Survey - Haffner, Reynolds, Madsen, et al. 2017, ApJS, in
COMMENT prep (as of 1/1/17).
COMMENT
COMMENT When citing these data, please state in your acknowledgements:
COMMENT
COMMENT "The Wisconsin H-Alpha Mapper and its Sky Survey have been funded
COMMENT primarily through awards from the U.S. National Science Foundation."
COMMENT
COMMENT Visit http://www.astro.wisc.edu/wham/ for updates, full documentation,
COMMENT and further detail about these data.
COMMENT
COMMENT To create this cube, the original one-degree-beam spectra are translated
COMMENT to a regularly-spaced velocity axis and then interpolated onto a
COMMENT regularly-spaced spatial grid using the Natural Neighbor method. See the
COMMENT website for the original, unmodified observations and other dataset
COMMENT formats.
COMMENT
COMMENT Data values are in units of milli-Rayleighs / (m/s) or, equivalently,
COMMENT Rayleighs / (km/s). 1 R = 10^6 / 4*pi photons cm^-2 s^-1 sr^-1 or 2.41 x
COMMENT 10^-7 erg cm^-2 s^-1 sr^-1 at H-Alpha. This intensity corresponds to an
COMMENT EM of 2.25 cm^-6 pc when T = 8000 K and assuming no extinction.
HISTORY PUTAST: Sep 12 13:28:52 2017 World Coordinate System parameters written
fn = '/home/justin/Research/WHAM/wham-ss-DR1-v161116-170912.fits'
whamdat1 = fits.getdata(fn)
whamhdr1 = fits.getheader(fn)
whamhdr1
SIMPLE = T / Written by IDL: Tue Sep 12 13:21:42 2017
BITPIX = 8 /
NAXIS = 0 /
EXTEND = T / File contains extensions
DATE = '2017-09-12' /
TAG = 'DR1 ' / General release tag
VERSION = '161116 ' / Generated from wham-ss-161116.sav
COMMENT WHAM Sky Survey - Release DR1 (161116)
COMMENT
COMMENT REFERENCES:
COMMENT
COMMENT Northern Sky Survey - Haffner, Reynolds, Tufte, Madsen, Jaehnig, &
COMMENT Percival. 2003, ApJS, 149, 405, doi:10.1086/378850.
COMMENT
COMMENT (Full) Sky Survey - Haffner, Reynolds, Madsen, et al. 2017, ApJS, in
COMMENT prep (as of 1/1/17).
COMMENT
COMMENT When citing these data, please state in your acknowledgements:
COMMENT
COMMENT "The Wisconsin H-Alpha Mapper and its Sky Survey have been funded
COMMENT primarily through awards from the U.S. National Science Foundation."
COMMENT
COMMENT Visit http://www.astro.wisc.edu/wham/ for updates, full documentation,
COMMENT and further detail about these data.
COMMENT
COMMENT Each row in this dataset represents an H-Alpha observation within a
COMMENT one-degree diameter beam centered at the location specified covering
COMMENT approximately -100 to +100 km/s (LSR) with 12 km/s spectral
COMMENT resolution. This dataset contains the original spectral observations
COMMENT and has not been resampled/regridded in location or velocity. Actual
COMMENT velocity coverage changes from row to row and observations are taken
COMMENT on an irregular grid. See the website for other dataset formats.
COMMENT
COMMENT Logitude (GAL-LON) and latitude (GAL-LAT) are given as degrees in the
COMMENT IAU (1958) Galactic reference frame.
COMMENT
COMMENT Velocities (VELOCITY) are given as km/s in the Local Standard of Rest
COMMENT reference frame using the conventional vector for peculiar solar
COMMENT motion of +20 km/s toward RA (1900) = 18h, DEC (1900) = +30d.
COMMENT
COMMENT Intensities (DATA), variances (VARIANCE), integrated errors (ERROR)
COMMENT are specified in Rayleighs (10^6 / 4*pi photons cm^-2 s^-1 sr^-1). 1 R
COMMENT = 2.41 x 10^-7 erg cm^-2 s^-1 sr^-1 at H-Alpha and is equivalent to an
COMMENT EM of 2.25 cm^-6 pc when T = 8000 K and assuming no extinction. The
COMMENT variances and errors listed here are only the intrinsic measurement
COMMENT error of the data points. They do not include the additional
COMMENT uncertainties associated with our absolute calibration or the
COMMENT systematic errors associated with the removal of the geocoronal and
COMMENT other atmospheric emission lines from the spectra. See the
COMMENT documentation at our website for more information.
COMMENT
COMMENT The integrated emission between -80 and +80 km/s (INTEN) and error
COMMENT (ERROR) are included for each row. These values are identical to the
COMMENT integrated dataset release (wham-ss-*-int.fits). If the ERROR in a row
COMMENT is zero, the original pointing is contaminated by a bright star (STAR
COMMENT = 1) with V(mag) < 6. The INTEN column represents the average of
COMMENT (good) nearby neighbor points. The original intensity (OINTEN) is
COMMENT provided for these pointings as a reference.
COMMENT
COMMENT
COMMENT BLOCK and POINTING numbers are an internal WHAM designation that
COMMENT groups and tracks observations. Our web page and publications may
COMMENT refer to these numbers when discussing systematic issues.
fn = '/home/justin/Research/WHAM/wham-ss-DR1-v161116-170912-grid.fits'
whamdat = fits.getdata(fn)
whamhdr = fits.getheader(fn)
Velocity range: -169.1 to 157.3 km/s
#masking for local velocities
velMask = (velAxis >= -80) & (velAxis <= 80)
velMask.sum()
np.int64(79)
whamIntegrated = np.nansum(whamdat[velMask, :, :] * abs(cdelt3), axis=0) / 1000
plt.imshow(whamIntegrated, origin = 'lower')
<matplotlib.image.AxesImage at 0x792897043050>
print(f"NaN count: {np.sum(np.isnan(whamdat))}")
print(f"Zero count: {np.sum(whamdat == 0)}")
print(f"Data range: {np.nanmin(whamdat):.4f} to {np.nanmax(whamdat):.4f}")
print(f"Integrated map range: {np.nanmin(whamIntegrated):.4f} to {np.nanmax(whamIntegrated):.4f}")
NaN count: 65383544
Zero count: 0
Data range: -0.0994 to 77.8963
Integrated map range: -1.3497 to 4417.2949
print(f"Channels selected: {velMask.sum()}")
print(f"Velocity range used: {velAxis[velMask][0]:.1f} to {velAxis[velMask][-1]:.1f} km/s")
Channels selected: 79
Velocity range used: -79.9 to 78.2 km/s
vmin = np.nanpercentile(whamIntegrated, 5)
vmax = np.nanpercentile(whamIntegrated, 95)
plt.figure(figsize=(12, 6))
plt.imshow(whamIntegrated, origin='lower', vmin=vmin, vmax=vmax, cmap='inferno')
plt.colorbar(label='H-alpha intensity (R)')
plt.title('WHAM Integrated H-alpha')
plt.show()
import matplotlib.pyplot as plt
import numpy as np
lon = np.linspace(180, -180, whamIntegrated.shape[1])
lat = np.linspace(-90, 90, whamIntegrated.shape[0])
lonGrid, latGrid = np.meshgrid(np.deg2rad(lon), np.deg2rad(lat))
vmin = np.nanpercentile(whamIntegrated, 5)
vmax = np.nanpercentile(whamIntegrated, 95)
fig = plt.figure(figsize=(12, 6))
ax = fig.add_subplot(111, projection='mollweide')
img = ax.pcolormesh(lonGrid, latGrid, whamIntegrated, cmap='inferno', vmin=vmin, vmax=vmax)
plt.colorbar(img, ax=ax, label='H-alpha Intensity (R)', shrink=0.5)
ax.set_xlabel('Galactic Longitude (deg)')
ax.set_ylabel('Galactic Latitude (deg)')
ax.grid(True, alpha=0.3)
plt.title('WHAM Integrated H-alpha Emission')
plt.tight_layout()
plt.show()
