In [2]:
# hack to get the path right
import sys
sys.path.append('..')
In [4]:
import ztf_sim
from astropy.time import Time
import pandas as pd
import numpy as np
import astropy.units as u
import pylab as plt
First we'll generate a test field grid. You only need to do this the first time you run the simulator.
In [7]:
ztf_sim.fields.generate_test_field_grid()
Let's load the Fields object with the default field grid. Fields is a thin wrapper around a pandas DataFrame containing the field information.
In [3]:
f = ztf_sim.fields.Fields()
The raw fieldid and coordinates are stored as a pandas Dataframe in the .fields attribute:
In [5]:
f.fields.head()
Out[5]:
Now let's calculate their altitude and azimuth at a specific time using the astropy.time.Time object:
In [7]:
f.alt_az(Time.now()).head()
Out[7]:
Demonstrating low-level access to fields by the fieldid index (usually not required):
In [5]:
f.fields.loc[853]
Out[5]:
We can select fields with conditionals:
In [8]:
f.fields['dec'] > -30.
Out[8]:
It's easier to use the select_fields convenience function, though. It returns a boolean Series indexed by fieldid that we can use to do calculations on subsets of the field grid.
In [15]:
cuts = f.select_fields(dec_range=[0,10],gridid=0,ecliptic_lat_range=[-5,5])
cuts.head()
Out[15]:
Calculate the current altitude and azimuth of the selected fields:
In [13]:
f.alt_az(Time.now(),cuts=cuts)
Out[13]:
Calculating the overhead time (max of ha, dec, dome slews and readout time):
In [12]:
f.overhead_time(853,Time.now())
Out[12]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
f = ztf_sim.fields.Fields()
Exposure_time = 60*u.second
Night_length=9*u.h
time0 = Time('2015-09-10 20:00:00') + 7*u.h
time = time0
f.fields = f.fields.join(pd.DataFrame(np.zeros(len(f.fields)),columns=['observed']))
f.fields = f.fields.join(pd.DataFrame(np.zeros(len(f.fields)),columns=['possibleToObserve']))
def observe(f, nightStart):
time=nightStart
goodAltitude = f.alt_az(time)['alt'] > 20
shouldObserve = f.fields['observed'] == 0
good = goodAltitude & shouldObserve #& f.alt_az(time+1*u.h)['alt'] < 20 # start with a field which won't be observable later
if np.all(good) is False:
good = goodAltitude & shouldObserve
fid = f.fields[good].iloc[0].name
f.fields['observed'][fid]+=1
f.fields['possibleToObserve'][goodAltitude] = 1
time += Exposure_time
while time < nightStart + Night_length:
goodAltitude = f.alt_az(time)['alt'] > 20
shouldObserve = f.fields['observed'] == 0
good = goodAltitude & shouldObserve
f.fields['possibleToObserve'][goodAltitude] = 1
if not np.any(good):
time += 60*u.s
continue
slewTime = f.overhead_time(fid,time)[good]
fid = int(slewTime.idxmin())
# print slewTime['overhead_time'][fid]
time += Exposure_time + slewTime['overhead_time'][fid]*u.second
f.fields['observed'][fid]+=1
# print time-7*u.h
# First night
observe(f,time)
fieldsPossible = np.sum(f.fields['possibleToObserve'])
print fieldsPossible
fieldsObserved = np.sum(f.fields['observed'])
print fieldsObserved
meanTime = (Night_length.to(u.s)-fieldsObserved*Exposure_time)/(fieldsObserved-1)
print meanTime
# Second night
time=time0+24*u.h
observe(f,time)
fieldsPossible = np.sum(f.fields['possibleToObserve'])
print fieldsPossible
fieldsObserved = np.sum(f.fields['observed'])
print fieldsObserved
meanTime = (2*Night_length.to(u.s)-fieldsObserved*Exposure_time)/(fieldsObserved-1)
print meanTime
In [ ]:
for dec in np.append(np.linspace(-90,90,10),0):
ra=np.linspace(0, 360,1000)
x,y = raDec2xy(ra,dec)
plt.plot(x,y,'k')
for ra in np.linspace(0,360,10):
dec=np.linspace(-90, 90,1000)
x,y = raDec2xy(ra,dec)
plt.plot(x,y,'k')
x,y = raDec2xy(f.fields['ra'],f.fields['dec'])
plt.plot(x,y,'o',color=(.8,.8,.8))
plt.show()
In [6]:
def raDec2xy(ra,dec):
# Using Aitoff projections (from Wiki) returns x-y coordinates on a plane of RA and Dec
theta = np.deg2rad(dec)
phi = np.deg2rad(ra)-np.pi #the range is [-pi,pi]
alpha=np.arccos(np.cos(theta)*np.cos(phi/2))
x=2*np.cos(theta)*np.sin(phi/2)/np.sinc(alpha/np.pi) # The python's sinc is normalazid, hence the /pi
y=np.sin(theta)/np.sinc(alpha/np.pi)
return x,y
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]: