Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Combined imaging/spectroscopy FoM #2

Open
wants to merge 3 commits into
base: main
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
113 changes: 98 additions & 15 deletions webserver/src/prism.py
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
import sys
import math
import numpy as np

def prism( inparams ):
# Survey
Expand All @@ -18,9 +19,12 @@ def prism( inparams ):
'tfixspec1': 900.0, # Exposure time of spectroscopic survey
'tfixspec2': 3600.0,
'tfixspec3': 0,
'z1': 2.0, # Redshift limit -- what does this mean? Complete?
'z1': 2.0, # Imaging redshift limit -- what does this mean? Complete?
'z2': 2.0,
'z3': 0.0,
'z1spec': 2.0, # Spectroscopy redshift limit
'z2spec': 2.0,
'z3spec': 0.0,
'nlim': 20, # Number of redshift bins from snflux ; don't change this!
# Instrument
'eff': 0.9, # Efficiency of... WHAT?
Expand All @@ -40,6 +44,7 @@ def prism( inparams ):
f5 = open('sninvarim.txt', 'w')
f4 = open('snfluxes.d', 'r')
f7 = open('errout312.txt', 'w')
f11 = open('sninvar.txt', 'w')

years=[params['length1'], params['length2'], params['length3']]

Expand All @@ -57,6 +62,10 @@ def prism( inparams ):
zmin=params['z1']
zmid= params['z2']
zmax=params['z3']
zminspec = params['z1spec']
zmidspec = params['z2spec']
zmaxspec = params['z3spec']

nlim=params['nlim']

eff=params['eff']
Expand Down Expand Up @@ -86,22 +95,24 @@ def prism( inparams ):

# Instrument parameters

pix=14.6
dc=0.015
pixim=8.97
dcim=0.015
zodiim=0.14
rnoise=5.0
res=100.0
resifu=80.0
wifu=3.0
lmin=1.0
lmax=2.0
rlamb=4.0
telarea=0.78*(math.pi*(120.0**2))
telarea=telarea*0.70
pixarea=0.11**2

# Old parameters, now unused
# pix = 14.6
# dcim = 0.015
# rnoise=5.0
# res=100.0
# resifu=80.0
# lmin=1.0
# lmax=2.0
# rlamb=4.0

rifu = [0,0,0,0,
580.0,
370.0,
Expand Down Expand Up @@ -189,7 +200,7 @@ def prism( inparams ):


snesp=[0.0]*nlim
nsnesp=[0.0]*nlim
# nsnesp=[0.0]*nlim
snmin=[0.0]*nlim
snmid=[0.0]*nlim
snmax=[0.0]*nlim
Expand All @@ -199,11 +210,11 @@ def prism( inparams ):

for i in range(1, nlim):
z=i*0.1+0.05
if z <= zmin:
if z <= zminspec:
snmin[i]=sr[i]*sqdgyrminsp*eff
if z <= zmid:
if z <= zmidspec:
snmid[i]=sr[i]*sqdgyrmidsp*eff
if z <= zmax:
if z <= zmaxspec:
snmax[i]=sr[i]*sqdgyrmaxsp*eff

sntot=snmin[i]+snmid[i]+snmax[i]
Expand Down Expand Up @@ -293,6 +304,7 @@ def prism( inparams ):
sigsys = [0.]*nlim
sigout = [0.]*nlim
sigdist = [0.]*nlim

for i in range(1, nlim):
z=i*0.1+0.05

Expand Down Expand Up @@ -360,7 +372,7 @@ def prism( inparams ):
z=i*0.1+0.05
nsne[i]=snesp[i]

sigstat[0]=0.17
sigstat[0]=0.13
sigsys[0]=0.0

f3.write(f'{om}\t{nw0}\t{nwa}\n')
Expand All @@ -369,6 +381,14 @@ def prism( inparams ):
f3.write('\t'.join([f'{x:.4f}' for x in sigstat])+'\n')
f3.write('\t'.join([f'{x:.4f}' for x in sigsys])+'\n')

# Store errors for combining surveys

nsnesp1=nsne.copy()
sigstatsp=[x/np.sqrt(y) for x, y in zip(sigstat,snesp)]
sigsyssp=sigsys.copy()
sigtotsp=[math.sqrt(x**2+y**2) for x, y in zip(sigstatsp, sigsyssp)]


# Imaging Survey. Xxxxxxxxxxxxxxxxxxxxxx

hrswide=((sqdegminim/0.28)*(tfixim[0]+70.0)*4.0)/(3600.0)
Expand Down Expand Up @@ -529,7 +549,7 @@ def prism( inparams ):
z=i*0.1+0.05
nsne[i]=sneim[i]

sigstat[0]=0.17
sigstat[0]=0.13
sigsys[0]=0.0

f5.write(f'{om}\t{nw0}\t{nwa}\n')
Expand All @@ -538,11 +558,74 @@ def prism( inparams ):
f5.write('\t'.join([f'{x:.4f}' for x in sigstat])+'\n')
f5.write('\t'.join([f'{x:.4f}' for x in sigsys])+'\n')

# Store errors for combining surveys

nsneim1=nsne.copy()
sigstatim=[x/np.sqrt(y) for x, y in zip(sigstat,sneim)]
sigsysim=sigsys.copy()
sigtotim=[math.sqrt(x**2 + y**2) for x, y in zip(sigstatim,sigsysim)]

# Combine spectroscopic and imaging surveys

f2.write(f'\nThe combined spectroscopic and imaging survey\n')
f2.write(f'\nz,nlow,nmid,nhigh,ntot,staterr,sigstat,sigsys,sigtotal\n')

nsnecomb = [x + y for x, y in zip(nsnesp1, nsneim1)]
snecomb = nsnecomb.copy()
errsum = [1/x**2 + 1/y**2 for x, y in zip(sigstatsp, sigstatim)]
sigstatcomb = 1/np.sqrt(errsum)
sigstatmulti = sigstatcomb*np.sqrt(snecomb)

# Calculate systematic errors a la Padmanabhan formulas

sigtotcomb = [0]*nlim
sigsyscomb = [0]*nlim
for i in range(nlim):

r=(sigsysim[i]*sigsyssp[i])/(sigtotim[i]*sigtotsp[i])
w1=(sigtotsp[i]**2-r*sigtotim[i]*sigtotsp[i])
w2=(sigtotim[i]**2+sigtotsp[i]**2-2.0*r*sigtotim[i]*sigtotsp[i])
w=w1/w2
s1=(w**2)*sigtotim[i]**2+((1.0-w)**2)*sigtotsp[i]**2
s2=2*w*(1.0-w)*r*sigtotim[i]*sigtotsp[i]
sq=s1+s2
sigtotcomb[i]=math.sqrt(sq)

sigsyscomb[i]=math.sqrt(sigtotcomb[i]**2-sigstatcomb[i]**2)


z=i*0.1+0.05

nmin=snim[0][i]+snsp[0][i]
nmid=snim[1][i]+snsp[1][i]
nmax=snim[2][i]+snsp[2][i]
nsnecomb[0]=800
ntotal=nsnecomb[i]


sigstatmulti[19]=0.235
sigsyscomb[15]=sigsysim[15]
sigsyscomb[16]=sigsysim[16]
sigsyscomb[17]=sigsysim[17]
sigsyscomb[18]=sigsysim[18]
sigsyscomb[19]=0.025


f2.write(f'{z:.2f}\t{int(nmin)}\t{int(nmid)}\t{int(nmax)}\t{int(ntotal)}\t{sigstatmulti[i]:.4f}\t{sigstatcomb[i]:.4f}\t{sigsyscomb[i]:.4f}\t{sigtotcomb[i]:.4f}\n')


f11.write(f'{om}\t{nw0}\t{nwa}\n')
f11.write('\t'.join([f'{int(x)}' for x in nsnecomb])+'\n')
f11.write(f'{nsamestat}\t{nsamesys}\t{nyescmb}\n')
f11.write('\t'.join([f'{x:.4f}' for x in sigstatmulti])+'\n')
f11.write('\t'.join([f'{x:.4f}' for x in sigsyscomb])+'\n')

f2.close()
f3.close()
f5.close()
f4.close()
f7.close()
f11.close()

except Exception as ex:
sys.stderr.write( f"Exception running prism; params is {params}\n" )
Expand Down
27 changes: 23 additions & 4 deletions webserver/src/snvar.py
Original file line number Diff line number Diff line change
Expand Up @@ -79,16 +79,23 @@ def ddwafn(x, u, v, om, q):
ddwafn=ewfull*(math.log(x)-1.+1./x)/(x**3.+q/om*ewfull)**1.5
return ddwafn

def snvar( inparams ):
def snvar(inparams, infile='sninvar.txt', mode='w', survey='combined'):

params = {} # All parameters are passed in via input file from prism but I'll keep this for convenience just in case
params.update( inparams )

try:
f1 = open('sninvar.txt', 'r')
f2 = open('snoutput.txt', 'w')
f1 = open(infile, 'r')
f2 = open('snoutput.txt', mode)
f3 = open('snerrout.txt', 'w')

if survey=='combined':
f2.write('========== FOM calculations for combined imaging/spectroscopic survey ==========\n\n')
elif survey=='imaging':
f2.write('\n========== FOM calculations for imaging survey only ==========\n\n')
elif survey=='spec':
f2.write('\n========== FOM calculations for spectroscopic survey only ==========\n\n')

# Input line 1: matter density Omega_m, w0, wa (fiducial 0.28, -1, 0)
om, w0, wa = [float(x) for x in f1.readline().split()]
# Line 2: number of SN per 0.1 bin in redshift, put 0 for empty bins, put z<0.1 SN in 1st bin
Expand Down Expand Up @@ -273,6 +280,12 @@ def snvar( inparams ):
f2.write('\nValues of w0, wa, reduced Fisher matrix F11, F12, F22\n')
f2.write('\t'.join([str(x) for x in [w0, wa, matsm[0, 0], matsm[0, 1], matsm[1, 1]]])+'\n')

if survey=='spec' and params['detail']==1:
f2.write('\n========== Detailed output ==========\n\n')
with open('prismoutSept15.txt', 'r') as f101:
f2.write(f101.read())


f1.close()
f2.close()
f3.close()
Expand All @@ -281,5 +294,11 @@ def snvar( inparams ):
sys.stderr.write( f"Exception running snvar; params is {params}\n" )
raise

def snvar_combined( inparams ):
snvar(inparams, infile='sninvar.txt', mode='w', survey='combined')
snvar(inparams, infile='sninvarim.txt', mode='a', survey='imaging')
snvar(inparams, infile='sninvarsp.txt', mode='a', survey='spec')


if __name__ == "__main__":
snvar( {} )
snvar_combined( {} )
52 changes: 28 additions & 24 deletions webserver/static/fom.js
Original file line number Diff line number Diff line change
Expand Up @@ -11,32 +11,36 @@ var fom = {};
fom.Context = function()
{

this.params = {'length1': { 'default': 1.9, 'desc': "Years of Survey 1" },
'length2': { 'default': 1.9, 'desc': "Years of Survey 2" },
'length3': { 'default': 1.9, 'desc': "Years of Survey 3" },
'sqdegim1': { 'default': 19.04, 'desc': "Square degrees of Survey 1" },
'sqdegim2': { 'default': 4.20, 'desc': "Square degrees of Survey 2" },
'sqdegim3': { 'default': 0.0, 'desc': "Square Degrees of Survey 3" },
'tfixim1': { 'default': 115.0, 'desc': "Imaging exposure time (s) of Survey 1" },
'tfixim2': { 'default': 450.0, 'desc': "Imaging exposure time (s) of Survey 2" },
'tfixim3': { 'default': 0, 'desc': "Imaging exposure time (s) of Survey 3" },
'z1': { 'default': 2.0, 'desc': "Redshift limit of Survey 1" },
'z2': { 'default': 2.0, 'desc': "Redshift limit of Survey 2" },
'z3': { 'default': 0.0, 'desc': "Redshift limit of Survey 3" },
this.params = {'length1': { 'default': 1.9, 'desc': "Years of Survey tier 1 (wide)" },
'length2': { 'default': 1.9, 'desc': "Years of Survey tier 2 (medium)" },
'length3': { 'default': 1.9, 'desc': "Years of Survey tier 3 (deep)" },
'sqdegim1': { 'default': 19.04, 'desc': "Square degrees of Survey tier 1 (wide)" },
'sqdegim2': { 'default': 4.20, 'desc': "Square degrees of Survey tier 2 (medium)" },
'sqdegim3': { 'default': 0.0, 'desc': "Square Degrees of Survey tier 3 (deep)" },
'tfixim1': { 'default': 115.0, 'desc': "Imaging exposure time (s) of Survey tier 1 (wide)" },
'tfixim2': { 'default': 450.0, 'desc': "Imaging exposure time (s) of Survey tier 2 (medium)" },
'tfixim3': { 'default': 0, 'desc': "Imaging exposure time (s) of Survey tier 3 (deep)" },
'z1': { 'default': 2.0, 'desc': "Imaging redshift limit of Survey tier 1 (wide)" },
'z2': { 'default': 2.0, 'desc': "Imaging redshift limit of Survey tier 2 (medium)" },
'z3': { 'default': 0.0, 'desc': "Imaging redshift limit of Survey tier 3 (deep)" },
// 'nlim': { 'default': 20, 'desc': "nlim(?)" },
'eff': { 'default': 0.9, 'desc': "Efficiency (frac of SNe caught)" },
'constsysim': { 'default': 0.015, 'desc': "Imaging systematic error at redshift divim" },
'divim': { 'default': 1.8, 'desc': "1+z where systematic error is constsysim" },
'stnim': { 'default': 10.0, 'desc': "S/N for imaging" },
'sqdegspec1': { 'default': 3.36, 'desc': "Square degrees for spectroscopy, Survey 1" },
'sqdegspec2': { 'default': 1.12, 'desc': "Square degrees for spectroscopy, Survey 2" },
'sqdegspec3': { 'default': 0.0, 'desc': "Square degrees for spectroscopy, Survey 3" },
'tfixspec1': { 'default': 900.0, 'desc': "Spectroscopy exposure time (s), Survey 1" },
'tfixspec2': { 'default': 3600.0, 'desc': "Spectroscopy exposure time (s), Survey 2" },
'tfixspec3': { 'default': 0, 'desc': "Spectroscopy exposure time (s), Survey 3" },
'constsysspec': { 'default': 0.01, 'desc': "Spectroscopy systematic error at redshit divspec" },
'divspec': { 'default': 1.8, 'desc': "1+z at which systematic error is constsysspec" },
'stnspec': { 'default': 20.0, 'desc': "Some sort of spectroscopic S/N..." },
'constsysim': { 'default': 0.015, 'desc': "Imaging systematic error coefficient" },
'divim': { 'default': 1.8, 'desc': "Imaging sysematic error denominator (1+z where systematic error is above)" },
'stnim': { 'default': 10.0, 'desc': "S/N for imaging to accept event" },
'sqdegspec1': { 'default': 3.36, 'desc': "Square degrees for spectroscopy, Survey tier 1 (wide)" },
'sqdegspec2': { 'default': 1.12, 'desc': "Square degrees for spectroscopy, Survey tier 2 (medium)" },
'sqdegspec3': { 'default': 0.0, 'desc': "Square degrees for spectroscopy, Survey tier 3 (deep)" },
'tfixspec1': { 'default': 900.0, 'desc': "Spectroscopy exposure time (s), Survey tier 1 (wide)" },
'tfixspec2': { 'default': 3600.0, 'desc': "Spectroscopy exposure time (s), Survey tier 2 (medium)" },
'tfixspec3': { 'default': 0, 'desc': "Spectroscopy exposure time (s), Survey tier 3 (deep)" },
'z1spec': { 'default': 2.0, 'desc': "Spectrscopy redshift limit of Survey tier 1 (wide)" },
'z2spec': { 'default': 2.0, 'desc': "Spectrscopy redshift limit of Survey tier 2 (medium)" },
'z3spec': { 'default': 0.0, 'desc': "Spectrscopy redshift limit of Survey tier 3 (deep)" },
'constsysspec': { 'default': 0.01, 'desc': "Spectroscopy systematic error coefficient" },
'divspec': { 'default': 1.8, 'desc': "spectroscopy sysematic error denominator (1+z where systematic error is above)" },
'stnspec': { 'default': 20.0, 'desc': "S/N for spectroscopy to accept event" },
'detail': { 'default': 0, 'desc': "Detailed output {0,1}" }
};
this.paramwidgets = {};
};
Expand Down
5 changes: 1 addition & 4 deletions webserver/webservice.py
Original file line number Diff line number Diff line change
Expand Up @@ -34,11 +34,8 @@ def calculate( logger=logging.getLogger("main") ):
logger.exception( f"Failed to run prism" )
return { 'status': 'error', 'error': f'Failed to run prism: {str(ex)}' }

sninvar = pathlib.Path( workdir / "sninvarim.txt" )
sninvar.replace( workdir / "sninvar.txt" )

try:
snvar.snvar( flask.request.json )
snvar.snvar_combined( flask.request.json )
except Exception as ex:
logger.exception( f"Failed to run snvar" )
return { 'status': 'error', 'error': f'Failed to run snvar: {str(ex)}' }
Expand Down