Wavelength recalibration with the sky lines (Mendel OH bands)

Introduction

Data can be recalibrated in wavelength easily when one OH sky lines are visible in most parts of the field of view. Cubes observed in the red part of the spectrum (e.g. SN3 filter) are especially interesting in this respect.

The general idea is to extract integrated spectra at different positions in the field of view and measure the velocity of the sky lines (which should be 0 if the calibration was perfect) (see Martin et al. 2017a).

skymap-fig0.svg

Then the correction map can be infered at each pixel by fitting a model.

skymap-fig1.svg

The calibration model is based on a simple modeling of the interferometer.

coords0.svg

First step: checking the calibration

[1]:
# import base class for the manipulation of a SITELLE spectral cube: HDFCube
from orcs.process import SpectralCube
import pylab as pl
[2]:
# load spectral cube
cube = SpectralCube('/home/thomas/M31_SN3.merged.cm1.1.0.hdf5')
master.fa09a|INFO| Cube is level 3
master.fa09a|INFO| shape: (2048, 2064, 840)
master.fa09a|INFO| wavenumber calibration: True
master.fa09a|INFO| flux calibration: True
master.fa09a|INFO| wcs calibration: True
[3]:
# extract and plot a spectrum of a large integrated region. Sky lines should be visible.
spectrum = cube.get_spectrum(50, 50, 20)
pl.figure(figsize=(10,6))
spectrum.plot()
pl.xlim(14500, 15500)
[3]:
(14500, 15500)
_images/script_example_wavelength_calibration_4_1.png

fitting the sky spectrum

we can also fit the integrated spectrum because there si a large number of sky lines and because we are not interested in a perfect fit (only the velocity is required), we can use a simple ‘sinc’ model with a fixed fwhm. All the lines are set to share the same velocity parameter and we set an inital guess of the velocity around 80 km/s. This general bias of 80 km/s is known (Martin et al. 2017b) and comes from the error made on the real wavelength of the calibration laser (which is falsely considered to be at 453.5 nm).

[4]:
sky_lines_cm1 = cube.get_sky_lines()
fit_res = spectrum.fit(sky_lines_cm1,
                       fmodel='sinc',
                       pos_def='1',
                       fwhm_def='fixed',
                       nofilter=False,
                       pos_cov=80)

print('Velocity: ', fit_res['velocity_gvar'][0])
pl.figure(figsize=(10,6))
spectrum.plot()
fit_res.get_spectrum().plot(c='red')
pl.xlim(14500, 15500)
Velocity:  77.5(1.5)
[4]:
(14500, 15500)
_images/script_example_wavelength_calibration_6_2.png

Second step: Mapping the sky velocity

This process, like other parallel processes may work better if run from a python script with the command. The following script can be used:

run_map_sky_velocity.py

from orcs.process import SpectralCube
if __name__ == "__main__": # this is very important when running parallel processes
    cube = SpectralCube('/home/thomas/M31_SN3.merged.cm1.1.0.hdf5')
    cube.map_sky_velocity(80, div_nb=10, exclude_reg_file_path='m31_exclude.reg', no_fit=False, plot=False) # the mean velocity bias is set around 80 km/s

Warning: this process can take a long time because a spectrum is extracted and fitted for each point of a 40x40 grid by default. For a basic example it is recommended to limit the number of grid division to 10x10.

[7]:
# once the fit is done, you can see the results with the following
cube = SpectralCube('/home/thomas/M31_SN3.merged.cm1.1.0.hdf5') # it is always better to reload the data cube right before running a parallel process
cube.map_sky_velocity(80, div_nb=10, exclude_reg_file_path='m31_exclude.reg', no_fit=False) # the mean velocity bias is set around 80 km/s
master.fa09a|INFO| Cube is level 3
master.fa09a|INFO| shape: (2048, 2064, 840)
master.fa09a|INFO| wavenumber calibration: True
master.fa09a|INFO| flux calibration: True
master.fa09a|INFO| wcs calibration: True
master.fa09a|WARNING| fitting process already done ! set no_fit=True if you do not want to redo it
master.fa09a|INFO| X range: 0 2048, Y range: 0 2064
master.fa09a|INFO| Radius: 30
master.fa09a|INFO| excluding region from file m31_exclude.reg
master.fa09a|INFO| 58 regions to fit
master.fa09a|INFO| 88 sky lines to fit
 loading region: Shape : circle ( Number(1861.81818),Number(1876.36364),Number(30.00000) )

master.fa09a|INFO| reloading cube
master.fa09a|INFO| Cube is level 3
master.fa09a|INFO| shape: (2048, 2064, 840)
master.fa09a|INFO| wavenumber calibration: True
master.fa09a|INFO| flux calibration: True
master.fa09a|INFO| wcs calibration: True
master.fa09a|INFO| Init of the parallel processing server with 32 threads
 [==========] [100%] [completed in 4m33s]
master.fa09a|INFO| Closing parallel processing server
master.fa09a|INFO| reloading cube

master.fa09a|INFO| Cube is level 3
master.fa09a|INFO| shape: (2048, 2064, 840)
master.fa09a|INFO| wavenumber calibration: True
master.fa09a|INFO| flux calibration: True
master.fa09a|INFO| wcs calibration: True
master.fa09a|INFO| parallel processing closed
master.fa09a|INFO| Velocity of the first line (km/s): 77.3(2.0)
master.fa09a|INFO| Velocity of the first line (km/s): 75.4(2.0)
master.fa09a|INFO| Velocity of the first line (km/s): 71.3(1.9)
master.fa09a|INFO| Velocity of the first line (km/s): 70.7(2.3)
master.fa09a|INFO| Velocity of the first line (km/s): 71.7(2.8)
master.fa09a|INFO| Velocity of the first line (km/s): 66.7(2.5)
master.fa09a|INFO| Velocity of the first line (km/s): 64.1(2.3)
master.fa09a|INFO| Velocity of the first line (km/s): 58.1(2.3)
master.fa09a|INFO| Velocity of the first line (km/s): 62.6(2.2)
master.fa09a|INFO| Velocity of the first line (km/s): 61.0(2.5)
master.fa09a|INFO| Velocity of the first line (km/s): 74.0(1.8)
master.fa09a|INFO| Velocity of the first line (km/s): 75.5(1.9)
master.fa09a|INFO| Velocity of the first line (km/s): 77.2(2.2)
master.fa09a|INFO| Velocity of the first line (km/s): 72.5(3.1)
master.fa09a|INFO| Velocity of the first line (km/s): 69.7(3.3)
master.fa09a|INFO| Velocity of the first line (km/s): 55.3(2.8)
master.fa09a|INFO| Velocity of the first line (km/s): 65.0(2.9)
master.fa09a|INFO| Velocity of the first line (km/s): 58.9(3.0)
master.fa09a|INFO| Velocity of the first line (km/s): 62.7(2.4)
master.fa09a|INFO| Velocity of the first line (km/s): 60.3(2.3)
master.fa09a|INFO| Velocity of the first line (km/s): 76.4(2.0)
master.fa09a|INFO| Velocity of the first line (km/s): 77.9(2.2)
master.fa09a|INFO| Velocity of the first line (km/s): 73.4(2.6)
master.fa09a|INFO| Velocity of the first line (km/s): 48.3(2.6)
master.fa09a|INFO| Velocity of the first line (km/s): 61.7(2.8)
master.fa09a|INFO| Velocity of the first line (km/s): 78.6(1.9)
master.fa09a|INFO| Velocity of the first line (km/s): 75.4(2.5)
master.fa09a|INFO| Velocity of the first line (km/s): 104.3(3.2)
master.fa09a|INFO| Velocity of the first line (km/s): 65.0(2.7)
master.fa09a|INFO| Velocity of the first line (km/s): 75.6(2.4)
master.fa09a|INFO| Velocity of the first line (km/s): 78.8(2.9)
master.fa09a|INFO| Velocity of the first line (km/s): 65.2(2.4)
master.fa09a|INFO| Velocity of the first line (km/s): 80.8(2.6)
master.fa09a|INFO| Velocity of the first line (km/s): 74.4(3.0)
master.fa09a|INFO| Velocity of the first line (km/s): 69.8(2.3)
master.fa09a|INFO| Velocity of the first line (km/s): 78.5(2.4)
master.fa09a|INFO| Velocity of the first line (km/s): 65.4(2.5)
master.fa09a|INFO| Velocity of the first line (km/s): 65.9(2.2)
master.fa09a|INFO| Velocity of the first line (km/s): 74.6(2.5)
master.fa09a|INFO| Velocity of the first line (km/s): 58.8(2.7)
master.fa09a|INFO| Velocity of the first line (km/s): 67.0(2.2)
master.fa09a|INFO| Velocity of the first line (km/s): 70.6(2.0)
master.fa09a|INFO| Velocity of the first line (km/s): 80.9(2.7)
master.fa09a|INFO| Velocity of the first line (km/s): 80.1(3.5)
master.fa09a|INFO| Velocity of the first line (km/s): 70.8(3.1)
master.fa09a|INFO| Velocity of the first line (km/s): 78.1(2.6)
master.fa09a|INFO| Velocity of the first line (km/s): 70.2(1.8)
master.fa09a|INFO| Velocity of the first line (km/s): 74.5(1.4)
master.fa09a|INFO| Velocity of the first line (km/s): 79.4(2.9)
master.fa09a|INFO| Velocity of the first line (km/s): 79.3(2.7)
master.fa09a|INFO| Velocity of the first line (km/s): 77.5(3.3)
master.fa09a|INFO| Velocity of the first line (km/s): 77.7(3.1)
master.fa09a|INFO| Velocity of the first line (km/s): 76.1(3.5)
master.fa09a|INFO| Velocity of the first line (km/s): 74.7(2.9)
master.fa09a|INFO| Velocity of the first line (km/s): 75.2(2.2)
master.fa09a|INFO| Velocity of the first line (km/s): 75.7(1.7)
master.fa09a|INFO| Velocity of the first line (km/s): 71.6(1.5)
master.fa09a|INFO| Velocity of the first line (km/s): 73.8(1.3)
master.fa09a|INFO| x: [ 186.18181818  186.18181818  186.18181818  186.18181818  186.18181818
  186.18181818  186.18181818  186.18181818  186.18181818  186.18181818
  372.36363636  372.36363636  372.36363636  372.36363636  372.36363636
  372.36363636  372.36363636  372.36363636  372.36363636  372.36363636
  558.54545455  558.54545455  558.54545455  558.54545455  558.54545455
  744.72727273  744.72727273  744.72727273  744.72727273  930.90909091
  930.90909091  930.90909091 1117.09090909 1117.09090909 1117.09090909
 1303.27272727 1303.27272727 1303.27272727 1489.45454545 1489.45454545
 1489.45454545 1489.45454545 1675.63636364 1675.63636364 1675.63636364
 1675.63636364 1675.63636364 1675.63636364 1861.81818182 1861.81818182
 1861.81818182 1861.81818182 1861.81818182 1861.81818182 1861.81818182
 1861.81818182 1861.81818182 1861.81818182]
master.fa09a|INFO| y: [ 187.63636364  375.27272727  562.90909091  750.54545455  938.18181818
 1125.81818182 1313.45454545 1501.09090909 1688.72727273 1876.36363636
  187.63636364  375.27272727  562.90909091  750.54545455  938.18181818
 1125.81818182 1313.45454545 1501.09090909 1688.72727273 1876.36363636
  187.63636364  375.27272727  562.90909091 1688.72727273 1876.36363636
  187.63636364  375.27272727 1688.72727273 1876.36363636  187.63636364
  375.27272727 1876.36363636  187.63636364 1688.72727273 1876.36363636
  187.63636364 1688.72727273 1876.36363636  187.63636364 1501.09090909
 1688.72727273 1876.36363636  187.63636364  375.27272727 1313.45454545
 1501.09090909 1688.72727273 1876.36363636  187.63636364  375.27272727
  562.90909091  750.54545455  938.18181818 1125.81818182 1313.45454545
 1501.09090909 1688.72727273 1876.36363636]
master.fa09a|INFO| v: [ 77.30690122  75.42723292  71.29603104  70.74703947  71.71676527
  66.74086248  64.11377386  58.053346    62.58297622  61.04474574
  74.01054084  75.52157459  77.22636344  72.5418529   69.69014005
  55.32754889  64.99891349  58.9375885   62.70848958  60.32799584
  76.36498222  77.90017749  73.41305969  48.25166759  61.66747683
  78.61262775  75.37893981 104.25308835  65.04902046  75.61380729
  78.84037233  65.16314474  80.77022818  74.38134783  69.81025738
  78.51941433  65.42526138  65.85118873  74.60799346  58.79111662
  66.95923754  70.63770242  80.9106735   80.12792716  70.77220517
  78.1282894   70.221003    74.52473849  79.4390245   79.2954783
  77.52699254  77.7411457   76.13299069  74.74027681  75.16744211
  75.65534916  71.59209525  73.79256937]
master.fa09a|WARNING| /home/thomas/Astro/Python/ORB/Orcs/orcs/utils.py:212: RuntimeWarning: invalid value encountered in greater
  vel[vel > np.nanpercentile(vel, 95)] = np.nan

master.fa09a|INFO| First laser wavelentgh calibration estimation: 543.3684614763755 nm
master.fa09a|INFO| Angle at the center of the frame: 15.494273258242744
master.fa09a|WARNING| /home/thomas/Astro/Python/ORB/Orb/orb/utils/image.py:1191: RuntimeWarning: invalid value encountered in less
  calib_laser_map[np.nonzero(calib_laser_map < value_min)] = np.nan

master.fa09a|INFO| > Binning calibration map
master.fa09a|WARNING| /home/thomas/Astro/Python/ORB/Orb/orb/utils/image.py:1074: RuntimeWarning: Mean of empty slice
  return np.squeeze(np.nanmean(np.nanmean(im_view, axis=3), axis=1))

master.fa09a|INFO| > Calibration laser map fit
master.fa09a|INFO|   > First fit on the central portion of the calibration laser map (30.0% of the total size)
master.fa09a|INFO|     > Calibration laser map fit parameters:
    distance to mirror: 23.681739781375743 cm
    X angle from the optical axis to the center: 0.0 degrees (Fixed)
    Y angle from the optical axis to the center: 15.501669386747308 degrees
    Tip-tilt angle of the detector along X: 0.0 degrees (Fixed)
    Tip-tilt angle of the detector along Y: 0.0 degrees (Fixed)
    Rotation angle of the detector: 0.0 degrees (Fixed)
    Calibration laser wavelength: 543.3684614763755 nm (Fixed)
    Error on fit: mean -1.2212972352744125e-06, std 0.049408481925003994 (in nm)
    Error on fit: mean -0.0006742923017409128, std 27.27899322170293 (in km/s)
master.fa09a|INFO|   > Second fit on the central portion of the calibration laser map (30.0% of the total size)
master.fa09a|INFO|     > Calibration laser map fit parameters:
    distance to mirror: 23.66483634860634 cm
    X angle from the optical axis to the center: -0.44808541110452904 degrees
    Y angle from the optical axis to the center: 15.495141935947101 degrees
    Tip-tilt angle of the detector along X: 0.41142114590561524 degrees
    Tip-tilt angle of the detector along Y: -1.41769820240921 degrees
    Rotation angle of the detector: 0.0 degrees (Fixed)
    Calibration laser wavelength: 543.3684614763755 nm (Fixed)
    Error on fit: mean -3.798378046820534e-07, std 0.0007481455641519504 (in nm)
    Error on fit: mean -0.0002097128366541575, std 0.41305980224865047 (in km/s)
master.fa09a|INFO|   > Third fit on a larger portion of the map (50.0% of the total size)
master.fa09a|INFO|     > Calibration laser map fit parameters:
    distance to mirror: 23.702051086393453 cm
    X angle from the optical axis to the center: -0.44853898242361845 degrees
    Y angle from the optical axis to the center: 15.495504360876263 degrees
    Tip-tilt angle of the detector along X: 0.25839602728725264 degrees
    Tip-tilt angle of the detector along Y: -0.8882041482610137 degrees
    Rotation angle of the detector: 0.0 degrees (Fixed)
    Calibration laser wavelength: 543.3684614763755 nm (Fixed)
    Error on fit: mean -2.4518726781407214e-06, std 0.003041168391745898 (in nm)
    Error on fit: mean -0.0013537072089970703, std 1.6790641750624251 (in km/s)
master.fa09a|INFO|   > Zernike polynomials fit of the residual wavefront
master.fa09a|WARNING| /home/thomas/Astro/Python/ORB/Orb/orb/ext/zern.py:300: FutureWarning: `rcond` parameter will change to the default of machine precision times ``max(M, N)`` where M and N are the input matrix dimensions.
To use the future default and silence this warning we advise to pass `rcond=None`, to keep using the old, explicitly pass `rcond=-1`.
  wf_zern_vec = np.linalg.lstsq((zern_basismat[:, grid_vec] * weight).T, wf_w.ravel())[0]

master.fa09a|INFO| Standard deviation of the residual: 1.8507972727041484e-05
master.fa09a|INFO| Median relative error (err/val)): 0.05%
master.fa09a|INFO| > final error (std on 50.0% of the total size): 1.763e-05 nm, 9.733e-03 km/s
master.fa09a|INFO|     > New calibration laser map fit parameters:
    distance to mirror: 23.67734145891187 cm
    X angle from the optical axis to the center: 0.1594504146076437 degrees
    Y angle from the optical axis to the center: 15.503103983197793 degrees
    Tip-tilt angle of the detector along X: 4.033689007793616 degrees
    Tip-tilt angle of the detector along Y: -0.3016934072811556 degrees
    Rotation angle of the detector: 2.2124798963045356 degrees
    Calibration laser wavelength: 543.3680110441915 nm

master.fa09a|INFO| fit residual std (in km/s):
master.fa09a|INFO| median error on the data (in km/s)
master.fa09a|INFO| Data written as ./M31_SN3/M31_SN3.calibration_laser_map.fits in 0.33 s
master.fa09a|INFO| Data written as ./M31_SN3/M31_SN3.wavefront_map.fits in 0.23 s
master.fa09a|INFO| Data written as ./M31_SN3/M31_SN3.skymap.fits in 0.26 s
_images/script_example_wavelength_calibration_8_7.png
_images/script_example_wavelength_calibration_8_8.png
_images/script_example_wavelength_calibration_8_9.png
_images/script_example_wavelength_calibration_8_10.png
_images/script_example_wavelength_calibration_8_11.png
_images/script_example_wavelength_calibration_8_12.png

Ouputs

The most important output is the skymap which gives the modelized velocity of the sky in the field of view and permit to correct the measured velocities in the cube.This map is named *skymap.fits. You can add it to the obtained velocity maps to correct the velocity.

In this particular case this is M31_SN3/M31_SN3.skymap.fits

The histogram above shows the difference between the measured velocity and the fitted velocity at each valid data point.

You can see that there are some data points at the center of the field which were excluded from the fit. This comes from the fact that M31 is very bright at the center and the sky lines are lost in the noise produced by the high continuum background. Modelling the interferometer makes possible to estimate the velocity of the lines at the center of the field where no velocity measurement could be done. And the result is stricking, it works particularly well even if the case of M31 (Martin et al. 2017a)

[14]:
import orb.utils.io
skymap = orb.utils.io.read_fits('./M31_SN3/M31_SN3.skymap.fits')
pl.imshow(skymap.T, origin='bottom')
pl.grid()
cb = pl.colorbar()
cb.set_label('velocity error (km/s)')
_images/script_example_wavelength_calibration_10_0.png
[ ]: