Interacting with Data Using SunPy Maps

In this example you will be learning how to create and modify SunPy Map objects.

Start by importing the necessary modules.

from __future__ import print_function, division

import numpy as np
import matplotlib.pyplot as plt

import astropy.units as u
from astropy.coordinates import SkyCoord

import sunpy.map
import sunpy.data.sample

SunPy Maps store 2D data in a numpy array and additional data in a metadata dictionary giving information relating to the data and instrument. You can create a Map in a number of ways, including loading a FITS file or URL: mymap = sunpy.map.Map('file1.fits') mymap = sunpy.map.Map(url_str) Or using creating manually by using tuple with the data/header within:

data = np.random.rand(20, 15)
header = {}
manual_map = sunpy.map.Map((data, header))

The data numpy array and metadata dictionary can easily be accessed:

print(manual_map.data)
print(manual_map.meta)

# In this case notice that the metadata has been populated by default with the
# naxis details that correspond to the array used for the data.

Out:

[[ 0.69836947  0.15828964  0.024308    0.21020753  0.33647155  0.97183899
   0.41852513  0.8118446   0.60884064  0.29032012  0.66228964  0.30482055
   0.00823923  0.00444755  0.74121314]
 [ 0.98417833  0.6463424   0.77323584  0.16218633  0.39769875  0.3451035
   0.75140441  0.42366572  0.38828395  0.47770309  0.95653026  0.59678239
   0.63100657  0.99008768  0.7963159 ]
 [ 0.77251345  0.1306636   0.28579398  0.16890246  0.61488043  0.80423626
   0.24804937  0.34361638  0.15247415  0.61078003  0.82854958  0.3201801
   0.54327189  0.0620441   0.39710567]
 [ 0.14771165  0.5334605   0.97535933  0.66023734  0.65533494  0.23110627
   0.82349009  0.21862137  0.77923525  0.09054379  0.45313036  0.86480797
   0.11572388  0.81295405  0.57542127]
 [ 0.09968168  0.39771611  0.39844202  0.0543224   0.77939374  0.81022548
   0.37440177  0.36546746  0.1647277   0.38983226  0.23649942  0.59749249
   0.90601487  0.53953533  0.59166391]
 [ 0.71487129  0.78271117  0.49376587  0.38366961  0.82082832  0.80297669
   0.28601195  0.66261399  0.42896906  0.66814888  0.44684566  0.82228401
   0.84015168  0.57334329  0.64658569]
 [ 0.76622077  0.02276971  0.80741366  0.98264962  0.21406299  0.54449126
   0.86398982  0.35186865  0.15632205  0.7708459   0.56087713  0.79307004
   0.34603063  0.79056122  0.58739867]
 [ 0.04121785  0.6398888   0.25506665  0.85013093  0.28709486  0.53755472
   0.06424225  0.39434257  0.65358809  0.21419187  0.90296039  0.49245644
   0.14244166  0.83242621  0.45125759]
 [ 0.11495204  0.50746005  0.3197024   0.47247831  0.49022621  0.16939177
   0.00826054  0.09382568  0.57405092  0.8290735   0.96703498  0.34680321
   0.75056808  0.54649843  0.58571733]
 [ 0.14935221  0.84208414  0.21446038  0.13961365  0.19122851  0.25364596
   0.66165724  0.67906991  0.65139115  0.99904835  0.76663766  0.73822859
   0.09118891  0.05582273  0.95515221]
 [ 0.98870044  0.76215308  0.03659258  0.19929276  0.94184655  0.01988983
   0.33598587  0.39613863  0.67552586  0.60597576  0.78763267  0.68242586
   0.16368761  0.02619077  0.49901267]
 [ 0.63257509  0.12816448  0.14340685  0.56707795  0.07540789  0.6751193
   0.58686566  0.08799346  0.0987314   0.28538872  0.81322899  0.24193317
   0.81278003  0.88511709  0.863526  ]
 [ 0.5657975   0.23438946  0.23535587  0.51255769  0.62519636  0.36396711
   0.82863932  0.52936452  0.72282017  0.75039641  0.7664685   0.26777324
   0.53882486  0.10402911  0.91059864]
 [ 0.7445537   0.3929592   0.26814312  0.6246745   0.5743438   0.39692059
   0.16741052  0.88278652  0.81498676  0.31259435  0.78829776  0.32602863
   0.42638347  0.3653615   0.02230023]
 [ 0.30192075  0.11541951  0.89171585  0.40140937  0.95325956  0.57928125
   0.41601544  0.33367804  0.86773127  0.59718125  0.99395878  0.35552307
   0.5154967   0.35951772  0.32434128]
 [ 0.87068832  0.40941669  0.05284311  0.91662836  0.74876191  0.28196544
   0.56961492  0.63991326  0.87814131  0.61729262  0.42289517  0.20023446
   0.08605663  0.03191522  0.81202126]
 [ 0.75258814  0.84805414  0.14078911  0.46208254  0.39130059  0.94524025
   0.53641506  0.69836329  0.58385512  0.97056303  0.27933935  0.24737607
   0.14351226  0.31168291  0.77681112]
 [ 0.14390819  0.0819527   0.74324414  0.52894979  0.12404246  0.76250024
   0.45856344  0.23049733  0.47461391  0.13007311  0.85502693  0.72835487
   0.56866337  0.10783609  0.84571397]
 [ 0.41985283  0.03143303  0.29757744  0.39072392  0.40336761  0.51351495
   0.90622741  0.13963322  0.54384428  0.35805351  0.94589852  0.83750067
   0.75354888  0.69827487  0.79211448]
 [ 0.25356988  0.94916745  0.64424828  0.02147514  0.88890316  0.99039698
   0.93534836  0.96209169  0.35874837  0.11257305  0.68909065  0.85628752
   0.46750214  0.47889936  0.99654006]]
MetaDict([('naxis1', 15), ('naxis2', 20), ('naxis', 2)])

You can quickly plot a map using the peek method:

manual_map.peek()
../../_images/sphx_glr_maps_example_001.png

SunPy Maps have a number of attributes that can be accessed easily, such as the x and y ranges:

print(manual_map.xrange)
print(manual_map.yrange)

# These return astropy Quantity objects.
# In general the attributes are populated using details in the metadata and in
# this case there is no centre pixel or pixel size information given so SunPy
# is defaulting to assuming each pixel is 1 arcsec.
# This is in Helioprojective tangent projection in both longitude and latitude:
print(manual_map.coordinate_system)

Out:

[ 0.49791667  0.50208333] deg
[ 0.49722222  0.50277778] deg
SpatialPair(axis1='HPLN-   ', axis2='HPLT-   ')

A real map example is given in the sample data, where the sunpy.data.sample.NAME returns the location of the given FITS file.

aia_map = sunpy.map.Map(sunpy.data.sample.AIA_171_IMAGE)
aia_map.peek(draw_limb=True)
../../_images/sphx_glr_maps_example_002.png

This has comprehensive metadata:

print(aia_map.meta)

Out:

MetaDict([('simple', True), ('bitpix', -32), ('naxis', 2), ('naxis1', 1024), ('naxis2', 1024), ('bld_vers', 'V5R12X'), ('lvl_num', 1.5), ('t_rec', '2011-06-07T06:33:03Z'), ('trecstep', 1.0), ('trecepoc', '1977.01.01_00:00:00_TAI'), ('trecroun', 1), ('origin', 'SDO'), ('date', '2012-10-16T19:45:34'), ('telescop', 'SDO'), ('instrume', 'AIA_3'), ('date-obs', '2011-06-07T06:33:02.77'), ('t_obs', '2011-06-07T06:33:02.88Z'), ('camera', 3), ('img_type', 'LIGHT'), ('exptime', 0.234256), ('expsdev', 0.000135), ('int_time', 0.507812), ('wavelnth', 171), ('waveunit', 'angstrom'), ('wave_str', '171_THIN'), ('fsn', 27194331), ('fid', 0), ('quallev0', 0), ('quality', 0), ('totvals', 16777216), ('datavals', 16777216), ('missvals', 0), ('percentd', 100.0), ('datamin', -8), ('datamax', 17722), ('datamedn', 18), ('datamean', 26.6888), ('datarms', 9429260.0), ('dataskew', 81.433), ('datakurt', 17076.7421875), ('datacent', 30.76), ('datap01', -1.0), ('datap10', 1.0), ('datap25', 3.0), ('datap75', 39.0), ('datap90', 64.0), ('datap95', 87.0), ('datap98', 127.0), ('datap99', 168.0), ('nsatpix', 29), ('oscnmean', 'nan'), ('oscnrms', 'nan'), ('flat_rec', 'aia.flatfield[:#30]'), ('nspikes', 977), ('ctype1', 'HPLN-TAN'), ('cunit1', 'arcsec'), ('crval1', 3.223099507700556), ('cdelt1', 2.402792), ('crpix1', 512.5), ('ctype2', 'HPLT-TAN'), ('cunit2', 'arcsec'), ('crval2', 1.385781353025793), ('cdelt2', 2.402792), ('crpix2', 512.5), ('crota2', -0.138829), ('r_sun', 1573.89688496), ('mpo_rec', 'sdo.master_pointing[:#394]'), ('inst_rot', 0.019327), ('imscl_mp', 0.599489), ('x0_mp', 2049.459961), ('y0_mp', 2049.030029), ('asd_rec', 'sdo.lev0_asd_0004[:#10672125]'), ('sat_y0', -4.519384), ('sat_z0', 13.793308), ('sat_rot', -3.4e-05), ('acs_mode', 'SCIENCE'), ('acs_eclp', 'NO'), ('acs_sunp', 'YES'), ('acs_safe', 'NO'), ('acs_cgt', 'GT3'), ('orb_rec', 'sdo.fds_orbit_vectors[2011.06.07_06:33:00_UTC]'), ('dsun_ref', 149597870691.0), ('dsun_obs', 151846026489.0), ('rsun_ref', 696000000.0), ('rsun_obs', 945.436711), ('gaex_obs', -13315441.72), ('gaey_obs', -25080995.82), ('gaez_obs', 31173320.91), ('haex_obs', -36634444009.2), ('haey_obs', -147360551017.0), ('haez_obs', 35731768.47), ('obs_vr', 86.972467), ('obs_vw', 31968.269521), ('obs_vn', 4842.822869), ('crln_obs', 340.698273), ('crlt_obs', 0.048591), ('car_rot', 2111), ('hgln_obs', 0.0), ('hglt_obs', 0.048591), ('roi_nwin', -2147483648), ('roi_sum', -2147483648), ('roi_nax1', -2147483648), ('roi_nay1', -2147483648), ('roi_llx1', -2147483648), ('roi_lly1', -2147483648), ('roi_nax2', -2147483648), ('roi_nay2', -2147483648), ('roi_llx2', -2147483648), ('roi_lly2', -2147483648), ('pixlunit', 'DN'), ('dn_gain', 17.7), ('eff_area', 2.419), ('eff_ar_v', 3.0), ('tempccd', -71.0), ('tempgt', 14.327), ('tempsmir', 34.029), ('tempfpad', 16.29), ('ispsname', 'aia.lev0_isp_0011'), ('isppktim', '2011-06-07T06:32:57.50Z'), ('isppktvn', '001.197'), ('aivnmst', 453), ('aimgots', 1686119616), ('asqhdr', 2174677979), ('asqtnum', 2), ('asqfsn', 27194331), ('aiahfsn', 27194323), ('aecdelay', 1537), ('aiaecti', 0), ('aiasen', 0), ('aifdbid', 241), ('aimgotss', 55109), ('aifcps', 10), ('aiftswth', 0), ('aifrmlid', 3057), ('aiftsid', 41729), ('aihismxb', 7), ('aihis192', 8377773), ('aihis348', 8384772), ('aihis604', 8386035), ('aihis860', 8386516), ('aifwen', 204), ('aimgshce', 237), ('aectype', 2), ('aecmode', 'ON'), ('aistate', 'CLOSED'), ('aiaecenf', 1), ('aifiltyp', 0), ('aimshobc', 41.104), ('aimshobe', 26.068001), ('aimshotc', 55.327999), ('aimshote', 69.344002), ('aimshcbc', 275.432007), ('aimshcbe', 260.484009), ('aimshctc', 289.556), ('aimshcte', 303.395996), ('aicfgdl1', 0), ('aicfgdl2', 107), ('aicfgdl3', 171), ('aicfgdl4', 236), ('aifoenfl', 1), ('aimgfsn', 5), ('aimgtyp', 0), ('aiawvlen', 7), ('aiagp1', 0), ('aiagp2', 0), ('aiagp3', 0), ('aiagp4', 0), ('aiagp5', 0), ('aiagp6', 0), ('aiagp7', 0), ('aiagp8', 619), ('aiagp9', 683), ('aiagp10', 748), ('agt1svy', 2), ('agt1svz', -6), ('agt2svy', 0), ('agt2svz', -1), ('agt3svy', -2), ('agt3svz', 3), ('agt4svy', -5), ('agt4svz', 5), ('aimgshen', 4), ('keywddoc', 'https:'), ('recnum', 76618218), ('blank', -32768), ('date_obs', '2011-06-07T06:33:02.77'), ('xcen', 2.91938326953), ('ycen', 1.09323792231), ('detector', 'AIA'), ('comment', "FITS (Flexible Image Transport System) format is defined in 'Astronomyand Astrophysics', volume 376, page 359; bibcode: 2001A&A...376..359HFITSHEAD2STRUCT"), ('history', 'FITSHEAD2STRUCT run at: Fri May 26 11:33:15 2017mreadfits_shm VERSION:1.20read_sdo VERSION:  2.10aia2wcsmin.pro VERSION:  5.10aia2wcsminaia2wcsmin  MPO_date: 2012-09-04T00:00:00Zaia2wcsmin  MPO_t_start: 2011-06-05T00:00:00Zaia2wcsmin  MPO_t_stop: 2011-06-12T00:00:00Zaia2wcsmin  MPO_version: 5ssw_reg.pro VERSION:  1.30ssw_regssw_reg  ROT called with cubic interpolation: cubic = -0.500000ssw_reg  Image registered to SDO image with FSN = 27194330ssw_reg  Image registered to SDO image with T_OBS = 2011-06-07T06:33:01.aia_fix_header.pro VERSION:  1.00aia_prep.pro VERSION: 5.10aia_reg.pro VERSION:  1.20'), ('keycomments', {'SIMPLE': 'conforms to FITS standard', 'NAXIS': 'number of array dimensions', 'BITPIX': 'array data type'})])

And find out information about the observation device and date:

print(aia_map.date)
print(aia_map.observatory)
print(aia_map.detector)
print(aia_map.exposure_time)
print(aia_map.coordinate_system)
print(aia_map.measurement)

Out:

2011-06-07 06:33:02.770000
SDO
AIA
0.234256 s
SpatialPair(axis1='HPLN-TAN', axis2='HPLT-TAN')
171.0 Angstrom

Maps also hold coordinate objects for the coordinate system they are in.

print(aia_map.coordinate_frame)

Out:

<Helioprojective Frame (obstime=2011-06-07 06:33:02.770000, rsun=696000000.0 m, observer=<HeliographicStonyhurst Coordinate (obstime=2011-06-07 06:33:02.770000): (lon, lat, radius) in (deg, deg, m)
    ( 0.,  0.048591,   1.51846026e+11)>)>

To see only a part of the image you create a submap, by specifying the top left and bottom right corners of the rectangle as either SkyCoord or Quantity objects.

bottom_left = aia_map.bottom_left_coord
top_right = SkyCoord(0 * u.arcsec, 0 * u.arcsec, frame=aia_map.coordinate_frame)
aia_submap = aia_map.submap(bottom_left, top_right)
aia_submap.peek(draw_limb=True)
../../_images/sphx_glr_maps_example_003.png

Similarly, if you want to reduce the angular resolution of the map you can use the resample method, specifying the dimensions as an Astropy Quantity in pixels:

dimensions = u.Quantity([50, 50], u.pixel)
aia_resampled_map = aia_map.resample(dimensions)
aia_resampled_map.peek(draw_limb=True, draw_grid=True)
../../_images/sphx_glr_maps_example_004.png

Similar to resampling you can use the superpixel method, this will reduce the resolution of the image by combining the number of pixels (in each dimension) in the dimensions argument into one single pixel. This can be used to increase the signal to noise ratio. For this the new dimensions must divide original image size exactly.

dimensions = u.Quantity(aia_map.dimensions) / 16
aia_superpixel_map = aia_map.superpixel(dimensions)
aia_superpixel_map.peek(draw_limb=True)
../../_images/sphx_glr_maps_example_005.png

Maps can also be rotated using the rotate method:

aia_rotated_submap = aia_submap.rotate(angle=10 * u.deg)
aia_rotated_submap.peek(draw_limb=True, draw_grid=True)

# Note: the data array is expanded so that none of the original data is lost
# through clipping.
plt.show()
../../_images/sphx_glr_maps_example_006.png

Total running time of the script: ( 0 minutes 3.117 seconds)

Generated by Sphinx-Gallery