notebook.community

Edit and run

GBM Geometry Demo

J. Michael Burgess

gbmeometry is a module with routines for handling GBM geometry. It performs a few tasks:

  • creates and astropy coordinate frame for Fermi GBM given a quarternion
  • allows for coordinate transforms from Fermi frame to in astropy frame (J2000, etc.)
  • plots the GBM NaI detectors at a given time for a given FOV
  • determines if a astropy SkyCoord location is within an NaI's FOV
  • creates interpolations over GBM quarternions and SC coordinates




In [1]:



    


%pylab inline
from astropy.coordinates import SkyCoord
import astropy.coordinates as coord
import astropy.units as u
from gbmgeometry import *



    


















    






Populating the interactive namespace from numpy and matplotlib

Making an interpolation from TRIGDAT

Getting the data

We can use GetGBMData to download the data





In [6]:



    


data = GetGBMData("080916009")


data.set_destination("") # You can enter a folder here. If you want the CWD, you do not have to set


data.get_trigdat()



    


















    






 [*********************100%***********************]  completed in 0.0 s

By default, GetGBMData will grad data from all detectors. However, one can set a subset to retrieve different data types





In [7]:



    


data.select_detectors('n1','n2','b0')


data.get_rsp_cspec()



    


















    






 [*********************100%***********************]  completed in 6.7 s

Interpolating

First let's create an interpolating object for a given TRIGDAT file (POSHIST files are also readable)





In [2]:



    


interp = PositionInterpolator(trigdat="glg_trigdat_all_bn080916009_v02.fit")



    











In [3]:



    


# In trigger times
print "Quaternions"
print interp.quaternion(0)
print interp.quaternion(10)
print
print "SC XYZ"
print interp.sc_pos(0)
print interp.sc_pos(10)



    


















    






Quaternions
[ 0.09894184  0.81399423  0.56763536  0.07357984]
[ 0.09651158  0.81315938  0.56970097  0.06998621]

SC XYZ
[ 3184.75  5985.5   1456.75]
[ 3111.77432458  6015.91372132  1488.98009345]

Single detector

One can look at a single detector which knows about it's orientation in the Fermi SC coordinates





In [4]:



    


na = NaIA(interp.quaternion(0))
print na.get_center()
print na.get_center().icrs #J2000
print na.get_center().galactic # Galactic
print 
print "Changing in time"
na.set_quaternion(interp.quaternion(100))
print na.get_center()
print na.get_center().icrs #J2000
print na.get_center().galactic # Galactic



    


















    






<SkyCoord (GBMFrame: sc_pos=None, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
    (123.73, -0.42)>
<SkyCoord (ICRS): (ra, dec) in deg
    (12.83264883, 51.9340734)>
<SkyCoord (Galactic): (l, b) in deg
    (122.91507726, -10.93767067)>

Changing in time
<SkyCoord (GBMFrame: sc_pos=None, quaternion=[ 0.07365115  0.80581046  0.58636083  0.03771755]): (Az, Zen) in deg
    (123.73, -0.42)>
<SkyCoord (ICRS): (ra, dec) in deg
    (14.14827804, 51.13038764)>
<SkyCoord (Galactic): (l, b) in deg
    (123.75793083, -11.73309177)>

We can also go back into the GBMFrame





In [5]:



    


center_j2000 = na.get_center().icrs
center_j2000



    


















    
Out[5]:








<SkyCoord (ICRS): (ra, dec) in deg
    (14.14827804, 51.13038764)>

















In [6]:



    


center_j2000.transform_to(GBMFrame(quaternion=interp.quaternion(100.)))



    


















    
Out[6]:








<SkyCoord (GBMFrame: sc_pos=None, quaternion=[ 0.07365115  0.80581046  0.58636083  0.03771755]): (Az, Zen) in deg
    (123.73020963, -0.4198826)>

Earth Centered Coordinates

The sc_pos are Earth centered coordinates (in km for trigdat and m for poshist) and can also be passed. It is a good idea to specify the units!





In [7]:



    


na = NaIA(interp.quaternion(0),interp.sc_pos(0)*u.km)
na.get_center()



    


















    
Out[7]:








<SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
    (123.73, -0.42)>

Working with the GBM class

Ideally, we want to know about many detectors. The GBM class performs operations on all detectors for ease of use. It also has plotting capabilities





In [8]:



    


myGBM = GBM(interp.quaternion(0),sc_pos=interp.sc_pos(0)*u.km)

myGBM.get_centers()



    


















    
Out[8]:








[<SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (123.73, -0.42)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (183.74, -0.32)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (236.61, 0.03)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (135.19, 44.45)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (0.0, 0.0)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (180.0, 0.0)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (45.89, 69.42)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (45.11, 44.69)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (58.44, -0.21)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (314.87, 44.76)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (303.15, -0.27)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (3.35, 0.03)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (224.93, 69.57)>,
 <SkyCoord (GBMFrame: sc_pos=[ 3184.75  5985.5   1456.75] km, quaternion=[ 0.09894184  0.81399423  0.56763536  0.07357984]): (Az, Zen) in deg
     (224.62, 43.82)>]

















In [9]:



    


[x.icrs for x in myGBM.get_centers()]



    


















    
Out[9]:








[<SkyCoord (ICRS): (ra, dec) in deg
     (12.83264883, 51.9340734)>, <SkyCoord (ICRS): (ra, dec) in deg
     (344.24965886, -2.97248228)>, <SkyCoord (ICRS): (ra, dec) in deg
     (318.5160438, -51.24281778)>, <SkyCoord (ICRS): (ra, dec) in deg
     (44.56221395, 13.56776347)>, <SkyCoord (ICRS): (ra, dec) in deg
     (165.84085686, -0.42750824)>, <SkyCoord (ICRS): (ra, dec) in deg
     (345.84085686, 0.42750824)>, <SkyCoord (ICRS): (ra, dec) in deg
     (90.02053912, -5.02944102)>, <SkyCoord (ICRS): (ra, dec) in deg
     (106.9639318, 13.09532427)>, <SkyCoord (ICRS): (ra, dec) in deg
     (137.09793576, 52.86180786)>, <SkyCoord (ICRS): (ra, dec) in deg
     (121.3144115, -45.95986919)>, <SkyCoord (ICRS): (ra, dec) in deg
     (194.28783066, -52.03023808)>, <SkyCoord (ICRS): (ra, dec) in deg
     (164.65699738, 2.70662314)>, <SkyCoord (ICRS): (ra, dec) in deg
     (58.29477735, -33.54731071)>, <SkyCoord (ICRS): (ra, dec) in deg
     (28.32615836, -45.28243674)>]

Plotting

We can look at the NaI view on the sky for a given FOV





In [10]:



    


myGBM.detector_plot(radius=60)



    


















    
























In [11]:



    


myGBM.detector_plot(radius=10,projection='ortho',lon_0=40)
myGBM.detector_plot(radius=10,projection='ortho',lon_0=0,lat_0=40,fignum=2)



    


















    






Warning: Cannot label parallels on Orthographic basemapWarning: Cannot label parallels on Orthographic basemap

In Fermi GBM coodinates

We can also plot in Fermi GBM spacecraft coordinates





In [12]:



    


myGBM.detector_plot(radius=60,fermi_frame=True)



    


















    
























In [13]:



    


myGBM.detector_plot(radius=10,projection='ortho',lon_0=20,fermi_frame=True)
myGBM.detector_plot(radius=10,projection='ortho',lon_0=200,fermi_frame=True,fignum=3)
myGBM.detector_plot(radius=10,projection='ortho',lon_0=0,lat_0=40,fignum=2,fermi_frame=True)



    


















    






Warning: Cannot label parallels on Orthographic basemapWarning: Cannot label parallels on Orthographic basemapWarning: Cannot label parallels on Orthographic basemap

Capturing points on the sky

We can even see which detector's FOVs contain a point on the sky. We create a mock GRB SKycoord first.





In [14]:



    


grb = SkyCoord(ra=130.,dec=-45 ,frame='icrs', unit='deg')


myGBM.detector_plot(radius=60,
                    projection='moll',
                    good=True, # only plot NaIs that see the GRB
                    point=grb,
                    lon_0=110,lat_0=-0)


myGBM.detector_plot(radius=60,
                    projection='ortho',
                    good=True, # only plot NaIs that see the GRB
                    point=grb,
                    lon_0=180,lat_0=-40,fignum=2)



    


















    






Warning: Cannot label parallels on Orthographic basemap

Looking at Earth occulted points on the sky

We can plot the points occulted by the Earth (assuming points 68.5 degrees form the Earth anti-zenith are hidden)





In [20]:



    


myGBM.detector_plot(radius=10,show_earth=True,lon_0=90)



    


















    
























In [16]:



    


myGBM.detector_plot(radius=10,lon_0=100,show_earth=True,projection='ortho')



    


















    






Warning: Cannot label parallels on Orthographic basemap









    
























In [17]:



    


myGBM.detector_plot(radius=10,show_earth=True,lon_0=120,lat_0=-30,fermi_frame=True,projection='ortho')



    


















    






Warning: Cannot label parallels on Orthographic basemap

Source/Detector Separation

We can even look at the separation angles for the detectors and the source





In [14]:



    


seps = myGBM.get_separation(grb)

seps.sort("Separation")

seps



    


















    
Out[14]:






<Table length=14>



DetectorSeparation


deg


str2float64


n36.16193940114


n441.7393116758


n053.0052176731


b054.6567413463


n654.8473106425


n556.7970862232


n161.7325232267


n766.3099675601


n883.475798479


n996.3850201623


n298.0510822102


nb123.163853651


b1125.343258654


na139.092756138

To see which detectors are valid, can look at the legal pairs map





In [13]:



    


get_legal_pairs()

Content source: drJfunk/gbmgeometry

Similar notebooks:

notebook.community | gallery | about