From 76d69a9236ab367f1d48755cdf4b5f42f45744f3 Mon Sep 17 00:00:00 2001
From: CI Bot <astrowisegitlab@rug.nl>
Date: Mon, 13 Feb 2023 12:31:19 +0000
Subject: [PATCH] Auto-format files with abstractautoformatter
---
notebooks/dark_raw_long.py | 86 ++++++++++
notebooks/dark_raw_short.py | 88 ++++++++++
notebooks/flat_raw_long.py | 124 ++++++++++++++
notebooks/flat_raw_short.py | 108 ++++++++++++
notebooks/illum_raw.py | 115 +++++++++++++
notebooks/sci_raw_and_calib.py | 218 ++++++++++++++++++++++++
notebooks/sci_stack.py | 84 +++++++++
notebooks/std_raw.py | 124 ++++++++++++++
src/micado/test/integration/cpl_test.py | 6 +-
9 files changed, 948 insertions(+), 5 deletions(-)
create mode 100644 notebooks/dark_raw_long.py
create mode 100644 notebooks/dark_raw_short.py
create mode 100644 notebooks/flat_raw_long.py
create mode 100644 notebooks/flat_raw_short.py
create mode 100644 notebooks/illum_raw.py
create mode 100644 notebooks/sci_raw_and_calib.py
create mode 100644 notebooks/sci_stack.py
create mode 100644 notebooks/std_raw.py
diff --git a/notebooks/dark_raw_long.py b/notebooks/dark_raw_long.py
new file mode 100644
index 000000000..d0d777a1e
--- /dev/null
+++ b/notebooks/dark_raw_long.py
@@ -0,0 +1,86 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: Python [conda env:scopesim]
+# language: python
+# name: conda-env-scopesim-py
+# ---
+
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+from astropy.time import Time, TimeDelta
+
+src = sim_tp.calibration.empty_sky()
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["IMG_4mas"],
+)
+
+# +
+start_time = datetime.datetime(2022, 1, 1, 21, 0, 0)
+
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 30,
+
+ "!OBS.catg": "CALIB",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "DARK",
+
+ "!OBS.filter_name_fw1": "close",
+ "!OBS.filter_name_fw2": "close",
+ "!OBS.filter_name_pupil": "close",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+ "!yamls": ["MICADO_H4RG.yaml"]
+})
+# -
+
+Time(start_time).to_value('mjd', 'long')
+
+for idark in range(3):
+ tm = TimeDelta(30.0 * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ micado = sim.OpticalTrain(cmd)
+
+ offeffect = micado.effects["element"] != "MICADO_DET"
+ for effect_name in list(micado.effects["name"][offeffect]):
+ try:
+ micado[effect_name].include = False
+ except:
+ #print (effect_name)
+ continue
+ for effect_name in ["filter_wheel_1", "filter_wheel_2", "pupil_wheel"]:
+ micado[effect_name].include = False
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+
+ micado.observe(src, seed=idark)
+ hdus = micado.readout()
+ hdus[0].writeto("dark_long_%d.fits"%idark, overwrite=True)
+
+hdus[0][0].header
+
+for effect_name in list(micado.effects["name"][micado.effects["included"]]):
+ print (effect_name, micado[effect_name].meta, "\n")
+
+
diff --git a/notebooks/dark_raw_short.py b/notebooks/dark_raw_short.py
new file mode 100644
index 000000000..4bac4fc18
--- /dev/null
+++ b/notebooks/dark_raw_short.py
@@ -0,0 +1,88 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: Python [conda env:scopesim]
+# language: python
+# name: conda-env-scopesim-py
+# ---
+
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+from astropy.time import Time, TimeDelta
+
+src = sim_tp.calibration.empty_sky()
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["IMG_4mas"],
+)
+
+# +
+start_time = datetime.datetime(2022, 1, 1, 21, 0, 0)
+
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 3,
+
+ "!OBS.catg": "CALIB",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "DARK",
+
+ "!OBS.filter_name_fw1": "close",
+ "!OBS.filter_name_fw2": "close",
+ "!OBS.filter_name_pupil": "close",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+ "!yamls": ["MICADO_H4RG.yaml"]
+})
+# -
+
+Time(start_time).to_value('mjd', 'long')
+
+for idark in range(3):
+ tm = TimeDelta(5.0 * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ micado = sim.OpticalTrain(cmd)
+
+ offeffect = micado.effects["element"] != "MICADO_DET"
+ for effect_name in list(micado.effects["name"][offeffect]):
+ try:
+ micado[effect_name].include = False
+ except:
+ #print (effect_name)
+ continue
+ for effect_name in ["filter_wheel_1", "filter_wheel_2", "pupil_wheel"]:
+ micado[effect_name].include = False
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+
+ micado.observe(src, seed=idark)
+ hdus = micado.readout()
+ hdus[0].writeto("dark_short_%d.fits"%idark, overwrite=True)
+
+hdus[0][0].header
+
+micado.effects
+
+for effect_name in list(micado.effects["name"][micado.effects["included"]]):
+ print (effect_name, micado[effect_name].meta, "\n")
+
+
diff --git a/notebooks/flat_raw_long.py b/notebooks/flat_raw_long.py
new file mode 100644
index 000000000..5d7196d11
--- /dev/null
+++ b/notebooks/flat_raw_long.py
@@ -0,0 +1,124 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: scopesim
+# language: python
+# name: scopesim
+# ---
+
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+from astropy.time import Time, TimeDelta
+
+src = sim_tp.calibration.flat_field()
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["IMG_4mas"],
+)
+
+# +
+start_time = datetime.datetime(2022, 1, 1, 22, 0, 0)
+
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 30,
+
+ "!OBS.catg": "CALIB",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "FLAT",
+
+ "!OBS.filter_name_fw1": "open",
+ "!OBS.filter_name_fw2": "H",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+ "!yamls": ["MICADO_H4RG.yaml"]
+})
+
+cmd["!DET.width"] = 4096 # pixel
+cmd["!DET.height"] = 4096
+# -
+
+cmd
+
+# +
+micado = sim.OpticalTrain(cmd)
+offeffect = micado.effects["element"] != "MICADO_DET"
+for effect_name in list(micado.effects["name"][offeffect]):
+ try:
+ micado[effect_name].include = False
+ except:
+ #print (effect_name)
+ continue
+
+micado["MICADO_DET"]["detector_window"].include = False
+micado["MICADO_DET"]["full_detector_array"].include = True
+# -
+
+micado.effects
+
+iflat = 0
+micado.observe(src, seed=iflat)
+hdus = micado.readout()
+hdus[0].writeto("flat_long_%d.fits"%iflat, overwrite=True)
+
+plt.figure(figsize=(20, 20))
+plt.imshow(hdus[0][1].data, origin="lower", norm=LogNorm())
+plt.colorbar()
+
+plt.figure(figsize=(20, 20))
+for hdu_i, plot_i in zip([1,2,3,6,5,4,7,8,9], range(1, 10)):
+ plt.subplot(3,3,plot_i)
+ plt.imshow(hdus[0][hdu_i].data, origin="lower", norm=LogNorm())
+
+hdus[0][0].header
+
+for iflat in range(1, 3):
+
+ tm = TimeDelta(30.0 * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ micado = sim.OpticalTrain(cmd)
+
+ offeffect = micado.effects["element"] != "MICADO_DET"
+ for effect_name in list(micado.effects["name"][offeffect]):
+ try:
+ micado[effect_name].include = False
+ except:
+ #print (effect_name)
+ continue
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+
+ micado.observe(src, seed=iflat)
+ hdus = micado.readout()
+ hdus[0].writeto("flat_long_%d.fits"%iflat, overwrite=True)
+
+for effect_name in list(micado.effects["name"][micado.effects["included"]]):
+ try:
+ print (effect_name, micado[effect_name].meta, "\n")
+ except:
+ continue
+
+for effect_name in ["filter_wheel_1", "filter_wheel_2", "pupil_wheel"]:
+ print (effect_name, micado[effect_name].meta, "\n")
+
+
+
diff --git a/notebooks/flat_raw_short.py b/notebooks/flat_raw_short.py
new file mode 100644
index 000000000..d9aaff48a
--- /dev/null
+++ b/notebooks/flat_raw_short.py
@@ -0,0 +1,108 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: scopesim
+# language: python
+# name: scopesim
+# ---
+
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+from astropy.time import Time, TimeDelta
+
+src = sim_tp.calibration.flat_field()
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["IMG_4mas"],
+)
+
+# +
+start_time = datetime.datetime(2022, 1, 1, 22, 0, 0)
+
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 3,
+
+ "!OBS.catg": "CALIB",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "FLAT",
+
+ "!OBS.filter_name_fw1": "open",
+ "!OBS.filter_name_fw2": "H",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+ "!yamls": ["MICADO_H4RG.yaml"]
+})
+
+cmd["!DET.width"] = 4096 # pixel
+cmd["!DET.height"] = 4096
+# -
+
+cmd
+
+# +
+micado = sim.OpticalTrain(cmd)
+offeffect = micado.effects["element"] != "MICADO_DET"
+for effect_name in list(micado.effects["name"][offeffect]):
+ try:
+ micado[effect_name].include = False
+ except:
+ #print (effect_name)
+ continue
+
+micado["MICADO_DET"]["detector_window"].include = False
+micado["MICADO_DET"]["full_detector_array"].include = True
+# -
+
+micado.effects
+
+for iflat in range(3):
+
+ tm = TimeDelta(3.0 * iflat * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ micado = sim.OpticalTrain(cmd)
+
+ offeffect = micado.effects["element"] != "MICADO_DET"
+ for effect_name in list(micado.effects["name"][offeffect]):
+ try:
+ micado[effect_name].include = False
+ except:
+ #print (effect_name)
+ continue
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+
+ micado.observe(src, random_seed=9001+iflat)
+ hdus = micado.readout()
+ hdus[0].writeto("flat_short_%d.fits"%iflat, overwrite=True)
+
+for effect_name in list(micado.effects["name"][micado.effects["included"]]):
+ try:
+ print (effect_name, micado[effect_name].meta, "\n")
+ except:
+ continue
+
+for effect_name in ["filter_wheel_1", "filter_wheel_2", "pupil_wheel"]:
+ print (effect_name, micado[effect_name].meta, "\n")
+
+
+
diff --git a/notebooks/illum_raw.py b/notebooks/illum_raw.py
new file mode 100644
index 000000000..20036f6a8
--- /dev/null
+++ b/notebooks/illum_raw.py
@@ -0,0 +1,115 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: scopesim
+# language: python
+# name: scopesim
+# ---
+
+# %load_ext autoreload
+# %autoreload 2
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+from astropy.time import Time, TimeDelta
+#from utils import fix_ext_header
+
+from scopesim_templates.micado import viking_fields as vf
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["IMG_4mas"],
+)
+
+# +
+start_time =datetime.datetime(2022, 1, 2, 1, 0, 0)
+shift_x = (0 * u.arcsec).to(u.degree)
+shift_y = (0 * u.arcsec).to(u.degree)
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 30,
+
+ "!OBS.catg": "CALIB",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "ILLUM",
+
+ "!OBS.filter_name_fw1": "open",
+ "!OBS.filter_name_fw2": "H",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+
+ # Fake dithering.
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": 1,
+
+})
+
+cmd["!DET.width"] = 4096 # pixel
+cmd["!DET.height"] = 4096
+
+# +
+xoffs = [3, -3, 3, -3]
+yoffs = [2, -2, -2, 2]
+
+for idither in range(4):
+
+ src = vf.viking_field(star_cat_id="illum", gal_cat_id="1", random_seed=9003)
+ # Shift the source in the opposite direction of the dither we want to simulate.
+ shift_x = (xoffs[idither] * u.arcsec).to(u.degree)
+ shift_y = (yoffs[idither] * u.arcsec).to(u.degree)
+ src.shift(-shift_x, -shift_y)
+
+ tm = TimeDelta(30.0 * idither * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ cmd.update(properties={
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": idither+1,
+ })
+ micado = sim.OpticalTrain(cmd)
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+ seeing_psf = sim.effects.psfs.SeeingPSF(fwhm = 1, name="psf")
+ micado.optics_manager.add_effect(seeing_psf)
+ micado["source_fits_keywords"].include = False
+ micado.observe(src, random_seed=idither)
+ hdus = micado.readout()
+ #fix_ext_header(hdus[0])
+ hdus[0].writeto("illum_%d.fits"%idither, overwrite=True)
+# -
+
+micado.effects
+
+hdus[0][0].header
+
+for effect_name in list(micado.effects["name"][micado.effects["included"]]):
+ try:
+ print (effect_name, micado[effect_name].meta, "\n")
+ except:
+ continue
+
+for effect_name in ["filter_wheel_1", "filter_wheel_2", "pupil_wheel"]:
+ print (effect_name, micado[effect_name].meta, "\n")
+
+
+micado.cmds.modes
+
+
diff --git a/notebooks/sci_raw_and_calib.py b/notebooks/sci_raw_and_calib.py
new file mode 100644
index 000000000..ad19ea82c
--- /dev/null
+++ b/notebooks/sci_raw_and_calib.py
@@ -0,0 +1,218 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: Python [conda env:scopesim]
+# language: python
+# name: conda-env-scopesim-py
+# ---
+
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+import numpy as np
+from astropy.wcs import WCS
+from astropy.io import fits
+from astropy.time import Time, TimeDelta
+from scopesim_templates.micado import viking_fields as vf
+#from utils import fix_ext_header
+
+shift_x = (0 * u.arcsec).to(u.degree)
+shift_y = (0 * u.arcsec).to(u.degree)
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["MCAO","IMG_4mas"],
+)
+
+# +
+start_time = datetime.datetime(2022, 1, 2, 3, 0, 0)
+
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 30,
+
+ "!OBS.catg": "SCIENCE",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "OBJECT",
+
+ "!OBS.filter_name_fw1": "open",
+ "!OBS.filter_name_fw2": "H",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+
+ # Fake dithering.
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": 1,
+
+})
+
+cmd["!DET.width"] = 4096 # pixel
+cmd["!DET.height"] = 4096
+# -
+
+# ## Raw Science
+
+# +
+micado = sim.OpticalTrain(cmd)
+
+micado["MICADO_DET"]["detector_window"].include = False
+micado["MICADO_DET"]["full_detector_array"].include = True
+# -
+
+micado.effects
+
+# +
+xoffs = [0, 3, -3, 3, -3]
+yoffs = [0, 2, -2, -2, 2]
+
+
+for idither in range(5):
+ src = vf.viking_field(star_cat_id="science", gal_cat_id="3", random_seed=9001)
+
+ # Shift the source in the opposite direction of the dither we want to simulate.
+ shift_x = (xoffs[idither] * u.arcsec).to(u.degree)
+ shift_y = (yoffs[idither] * u.arcsec).to(u.degree)
+ src.shift(-shift_x, -shift_y)
+
+ tm = TimeDelta(30.0 * idither * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ cmd.update(properties={
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": idither + 1,
+ })
+ micado = sim.OpticalTrain(cmd)
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+ micado["source_fits_keywords"].include = False
+ micado.observe(src, random_seed=9001+idither)
+ hdus = micado.readout()
+ #fix_ext_header(hdus[0])
+ hdus[0].writeto("sci_%d.fits"%idither, overwrite=True)
+# -
+
+# ## Calibrated
+
+# +
+# Turn off some effects for calibrated data
+
+micado_cal = sim.OpticalTrain(cmd)
+
+micado_cal["MICADO_DET"]["detector_window"].include = False
+micado_cal["MICADO_DET"]["full_detector_array"].include = True
+
+for effect_name in ["skycalc_atmosphere", "qe_curve",
+ "dark_current", "detector_linearity",
+ "readout_noise"]:
+ micado_cal[effect_name].include = False
+
+micado_sci = sim.OpticalTrain(cmd)
+
+micado_sci["MICADO_DET"]["detector_window"].include = False
+micado_sci["MICADO_DET"]["full_detector_array"].include = True
+
+# Turn off noise to calculate noise layer
+for effect_name in ["skycalc_atmosphere", "qe_curve",
+ "dark_current", "detector_linearity",
+ "readout_noise", "shot_noise"]:
+ micado_sci[effect_name].include = False
+# -
+
+micado_cal.effects
+
+micado_sci.effects
+
+# +
+xoffs = [0, 3, -3, 3, -3]
+yoffs = [0, 2, -2, -2, 2]
+
+for idither in range(5):
+ src = vf.viking_field(star_cat_id="science", gal_cat_id="3", random_seed=9001)
+
+ # Shift the source in the opposite direction of the dither we want to simulate.
+ shift_x = (xoffs[idither] * u.arcsec).to(u.degree)
+ shift_y = (yoffs[idither] * u.arcsec).to(u.degree)
+ src.shift(-shift_x, -shift_y)
+
+ tm = TimeDelta(30.0 * idither * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ cmd.update(properties={
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": idither + 1,
+ })
+
+
+ ## Set up optical trains
+ # Turn off some effects for calibrated data
+
+ micado_cal = sim.OpticalTrain(cmd)
+
+ micado_cal["MICADO_DET"]["detector_window"].include = False
+ micado_cal["MICADO_DET"]["full_detector_array"].include = True
+
+ for effect_name in ["skycalc_atmosphere", "qe_curve",
+ "dark_current", "detector_linearity",
+ "readout_noise"]:
+ micado_cal[effect_name].include = False
+
+ micado_sci = sim.OpticalTrain(cmd)
+
+ micado_sci["MICADO_DET"]["detector_window"].include = False
+ micado_sci["MICADO_DET"]["full_detector_array"].include = True
+ # Turn off noise to calculate noise layer
+ for effect_name in ["skycalc_atmosphere", "qe_curve",
+ "dark_current", "detector_linearity",
+ "readout_noise", "shot_noise"]:
+ micado_sci[effect_name].include = False
+
+
+ micado_cal.observe(src, random_seed=idither)
+ hdus = micado_cal.readout()
+
+ hdus[0][0].header.set("HIERARCH ESO PRO CATG", "IMG_CALIB")
+ #fix_hdus = fix_ext_header(hdus[0], calib = True)
+
+
+ # Add ERR and QC layers
+ micado_sci.observe(src, random_seed=idither)
+ hdus_pure = micado_sci.readout()
+
+ for hdu_i in range(1, 10):
+ hdus[0].insert(2 + 3*(hdu_i-1),
+ fits.ImageHDU(hdus[0][1 + 3*(hdu_i-1)].data-hdus_pure[0][hdu_i].data,
+ header = hdus[0][1 + 3*(hdu_i-1)].header))
+ hdus[0].insert(3 + 3*(hdu_i-1),
+ fits.ImageHDU((hdus[0][1 + 3*(hdu_i-1)].data-hdus[0][1 + 3*(hdu_i-1)].data).astype(int),
+ header = hdus[0][1 + 3*(hdu_i-1)].header))
+ for hdu_i in range(1, 10):
+ hdus[0][(hdu_i-1)*3 + 1].name = "DET%d.DATA"%hdu_i
+ hdus[0][(hdu_i-1)*3 + 2].name = "DET%d.ERR"%hdu_i
+ hdus[0][(hdu_i-1)*3 + 3].name = "DET%d.DQ"%hdu_i
+
+ hdus[0].writeto("sci_calib_%d.fits"%idither, overwrite=True)
+# -
+
+
diff --git a/notebooks/sci_stack.py b/notebooks/sci_stack.py
new file mode 100644
index 000000000..6c58f4a76
--- /dev/null
+++ b/notebooks/sci_stack.py
@@ -0,0 +1,84 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: Python [conda env:scopesim]
+# language: python
+# name: conda-env-scopesim-py
+# ---
+
+# +
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+
+import astropy.units as u
+import numpy as np
+from astropy.wcs import WCS
+from astropy.io import fits
+# -
+
+from reproject import reproject_exact, reproject_interp
+from reproject.mosaicking import reproject_and_coadd
+from reproject.mosaicking import find_optimal_celestial_wcs
+
+hdul = []
+img_hdul = []
+err_hdrl = []
+for idther in range(5):
+ sci_calib = fits.open("sci_calib_%d.fits"%idther)
+ hdul.append(sci_calib)
+ img_hdul += sci_calib[1::3]
+ err_hdrl += sci_calib[2::3]
+
+wcs_out, shape_out = find_optimal_celestial_wcs(img_hdul, resolution=0.004 * u.arcsec)
+
+array, footprint = reproject_and_coadd(img_hdul,wcs_out, shape_out=shape_out,
+ reproject_function= reproject_interp)
+
+array_err, __ = reproject_and_coadd(err_hdrl,wcs_out, shape_out=shape_out,
+ reproject_function= reproject_interp)
+
+np.nanmax(array), np.nanmean(array)
+
+plt.figure(figsize=(20, 20))
+ax1 = plt.subplot(projection=wcs_out)
+im1 = ax1.imshow(array, origin='lower', norm=LogNorm(vmax=100, vmin=10))
+ax1.set_title('Mosaic')
+
+plt.figure(figsize=(20, 20))
+ax1 = plt.subplot(projection=wcs_out)
+im1 = ax1.imshow(array_err, norm=LogNorm(vmax=20, vmin=2))
+ax1.set_title('Mosaic Error')
+
+plt.figure(figsize=(20, 20))
+ax2 = plt.subplot(projection=wcs_out)
+im2 = ax2.imshow(footprint, origin='lower')
+ax2.set_title('Footprint')
+
+st_hdu = hdul[0][0]
+
+st_hdul = fits.HDUList([st_hdu])
+st_hdul.append(fits.ImageHDU(array, header = st_hdul[0].header))
+st_hdul[0].header = hdul[0][0].header
+wcs_hd = wcs_out.to_header()
+for tt in wcs_hd:
+ st_hdul[1].header.set(tt, wcs_hd[tt])
+st_hdul[0].header.set("HIERARCH ESO PRO CATG", "IMG_STACK")
+
+st_hdul.append(fits.ImageHDU(array_err, header = st_hdul[1].header))
+st_hdul.append(fits.ImageHDU(footprint, header = st_hdul[1].header))
+st_hdul[1].name = "DATA"
+st_hdul[2].name = "ERR"
+st_hdul[3].name = "DQ"
+
+st_hdul.writeto('sci_stacked.fits', overwrite=True)
+
+
diff --git a/notebooks/std_raw.py b/notebooks/std_raw.py
new file mode 100644
index 000000000..3c5930a95
--- /dev/null
+++ b/notebooks/std_raw.py
@@ -0,0 +1,124 @@
+# ---
+# jupyter:
+# jupytext:
+# text_representation:
+# extension: .py
+# format_name: light
+# format_version: '1.5'
+# jupytext_version: 1.13.1
+# kernelspec:
+# display_name: Python [conda env:scopesim]
+# language: python
+# name: conda-env-scopesim-py
+# ---
+
+# %load_ext autoreload
+# %autoreload 2
+# %matplotlib inline
+from matplotlib import pyplot as plt
+from matplotlib.colors import LogNorm
+import datetime
+from astropy import units as u
+import scopesim as sim
+import scopesim_templates as sim_tp
+import astropy.units as u
+from astropy.time import Time, TimeDelta
+#from utils import fix_ext_header
+
+from scopesim_templates.micado import viking_fields as vf
+
+shift_x = (0 * u.arcsec).to(u.degree)
+shift_y = (0 * u.arcsec).to(u.degree)
+
+cmd = sim.UserCommands(
+ use_instrument="MICADO",
+ set_modes=["MCAO","IMG_4mas"],
+)
+
+# +
+start_time =datetime.datetime(2022, 1, 2, 2, 0, 0)
+
+cmd.update(properties={
+ "!OBS.ndit": 1,
+ "!OBS.dit": 30,
+
+ "!OBS.catg": "CALIB",
+ "!OBS.tech": "IMAGE",
+ "!OBS.type": "PHOTOM",
+
+ "!OBS.filter_name_fw1": "open",
+ "!OBS.filter_name_fw2": "H",
+
+ "!OBS.mjdobs": Time(start_time).mjd,
+ "!OBS.tplstart": start_time,
+
+ # Fake dithering.
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": 1,
+
+})
+
+cmd["!DET.width"] = 4096 # pixel
+cmd["!DET.height"] = 4096
+
+# +
+micado = sim.OpticalTrain(cmd)
+
+micado["MICADO_DET"]["detector_window"].include = False
+micado["MICADO_DET"]["full_detector_array"].include = True
+# -
+
+micado.effects
+
+# +
+xoffs = [3, -3, 3, -3]
+yoffs = [2, -2, -2, 2]
+
+for idither in range(4):
+
+ src = vf.viking_field(star_cat_id="stdstar", gal_cat_id="2", random_seed=9002)
+ # Shift the source in the opposite direction of the dither we want to simulate.
+ shift_x = (xoffs[idither] * u.arcsec).to(u.degree)
+ shift_y = (yoffs[idither] * u.arcsec).to(u.degree)
+ src.shift(-shift_x, -shift_y)
+
+ tm = TimeDelta(30.0 * idither * u.s)
+
+ cmd.update(properties={
+ "!OBS.mjdobs": (Time(start_time)+tm).mjd,
+ })
+
+ cmd.update(properties={
+ "!OBS.ra": shift_x,
+ "!OBS.dec": shift_y,
+ "!OBS.tplexpno": idither,
+ })
+ micado = sim.OpticalTrain(cmd)
+
+ micado["MICADO_DET"]["detector_window"].include = False
+ micado["MICADO_DET"]["full_detector_array"].include = True
+ micado["source_fits_keywords"].include = False
+ micado.observe(src, random_seed=9001+idither)
+ hdus = micado.readout()
+ #fix_ext_header(hdus[0])
+ hdus[0].writeto("std_%d.fits"%idither, overwrite=True)
+# -
+
+hdus[0][0].header
+
+for effect_name in list(micado.effects["name"][micado.effects["included"]]):
+ try:
+ print (effect_name, micado[effect_name].meta, "\n")
+ except:
+ continue
+
+for effect_name in ["filter_wheel_1", "filter_wheel_2", "pupil_wheel"]:
+ print (effect_name, micado[effect_name].meta, "\n")
+
+
+micado.cmds.modes
+
+micado["maory_generic_psf"]
+
+
diff --git a/src/micado/test/integration/cpl_test.py b/src/micado/test/integration/cpl_test.py
index f67259fb3..a0a89ff18 100644
--- a/src/micado/test/integration/cpl_test.py
+++ b/src/micado/test/integration/cpl_test.py
@@ -24,11 +24,7 @@ def test_cpl_processing():
"""Test whether CPL and prototype processing is identical."""
# TODO: do something smarter than marking the test as xfail, e.g.
# only run the test when the MICADO Pipeline is installed.
- (
- mymasterdark_proto,
- dark_real,
- mydarks_raw,
- ) = create_masterdark(
+ (mymasterdark_proto, dark_real, mydarks_raw,) = create_masterdark(
make=True,
ndetectors=9,
)
--
GitLab