Dithering

specsim/fiberloss.py

@@ -33,13 +33,17 @@ class GalsimFiberlossCalculator(object):
         Beta parameter value for the atmospheric PSF Moffat profile.
     maximum_fft_size : int
         Maximum size of FFT allowed.
+    fiber_placement: array
+        Fiber placement position in microns. It is needed for dithering
     """
     def __init__(self, fiber_diameter, wlen_grid, num_pixels=16,
-                 oversampling=32, moffat_beta=3.5, maximum_fft_size=32767):
+                 oversampling=32, moffat_beta=3.5, maximum_fft_size=32767,
+                 fiber_placement=[0., 0.]):
 
         self.wlen_grid = np.asarray(wlen_grid)
         self.moffat_beta = moffat_beta
-
+        self.fiber_placement = fiber_placement
+
         # Defer import to runtime.
         import galsim
 
@@ -48,21 +52,24 @@ def __init__(self, fiber_diameter, wlen_grid, num_pixels=16,
         # rather than on-sky angles.
         scale = fiber_diameter / num_pixels
         self.image = galsim.Image(num_pixels, num_pixels, scale=scale)
-
         self.gsparams = galsim.GSParams(maximum_fft_size=32767)
 
         # Prepare an anti-aliased image of the fiber aperture.
         nos = num_pixels * oversampling
-        dxy = (np.arange(nos) + 0.5 - 0.5 * nos) / (0.5 * nos)
-        rsq = dxy ** 2 + dxy[:, np.newaxis] ** 2
+
+        #dxy = (np.arange(nos) + 0.5 - 0.5 * nos) / (0.5 * nos)
+        dx = (np.arange(nos) + 0.5 - 0.5 * nos) / (0.5 * nos) - self.fiber_placement[0]/0.5/fiber_diameter
+        dy = (np.arange(nos) + 0.5 - 0.5 * nos) / (0.5 * nos) - self.fiber_placement[1]/0.5/fiber_diameter
+
+        #rsq = dxy ** 2 + dxy[:, np.newaxis] ** 2
+        rsq = dx ** 2 + dy[:, np.newaxis] ** 2
         inside = (rsq <= 1).astype(float)
         s0, s1 = inside.strides
         blocks = numpy.lib.stride_tricks.as_strided(
             inside, shape=(num_pixels, num_pixels, oversampling, oversampling),
             strides=(oversampling * s0, oversampling * s1, s0, s1))
         self.aperture = blocks.sum(axis=(2, 3)) / oversampling ** 2
 
-
     def create_source(self, fractions, half_light_radius,
                       minor_major_axis_ratio, position_angle):
         """Create a model for the on-sky profile of a single source.
@@ -241,7 +248,7 @@ def calculate(self, seeing_fwhm, scale, offset, blur_rms,
                 convolved.drawImage(offset=offsets, **draw_args)
                 # Calculate the fiberloss fraction for this fiber and wlen.
                 fiberloss[j, i] = np.sum(self.image.array * self.aperture)
-
+                
                 if saved_images_file is not None:
                     header['FIBER'] = j
                     header['WLEN'] = wlen
@@ -260,6 +267,7 @@ def calculate(self, seeing_fwhm, scale, offset, blur_rms,
                     header['COMMENT'] = 'Atmospheric seeing model'
                     hdu_list.append(astropy.io.fits.ImageHDU(
                         data=self.image.array.copy(), header=header))
+
                     if wlen == self.wlen_grid[-1]:
                         # Render the source profile without any offset after
                         # all other postage stamps for this fiber.
@@ -280,7 +288,8 @@ def calculate_fiber_acceptance_fraction(
     focal_x, focal_y, wavelength, source, atmosphere, instrument,
     source_types=None, source_fraction=None, source_half_light_radius=None,
     source_minor_major_axis_ratio=None, source_position_angle=None,
-    oversampling=32, saved_images_file=None, saved_table_file=None):
+    oversampling=32, saved_images_file=None, saved_table_file=None,
+    fiber_placement=[0., 0.]):
     """Calculate the acceptance fraction for a single fiber.
 
     The behavior of this function is customized by the instrument.fiberloss
@@ -345,7 +354,9 @@ def calculate_fiber_acceptance_fraction(
         extension determines the file format, and .ecsv is recommended.
         The saved file can then be used as a pre-tabulated input with
         instrument.fiberloss.method = 'table'.
-
+    fiber_placement: array
+        Fiber placement position in microns. It is needed for dithering
+
     Returns
     -------
     numpy array
@@ -384,7 +395,8 @@ def calculate_fiber_acceptance_fraction(
         wlen_grid.to(u.Angstrom).value,
         instrument.fiberloss_num_pixels,
         oversampling,
-        atmosphere.seeing_moffat_beta)
+        atmosphere.seeing_moffat_beta,
+        fiber_placement=fiber_placement)
 
     # Calculate the focal-plane optics at the fiber locations.
     scale, blur, offset = instrument.get_focal_plane_optics(

specsim/quickfiberloss.py

@@ -63,7 +63,7 @@ def main(args=None):
                             args.num_wlen) * wlen_unit
 
     # Calculate the seeing at each wavelength.
-    simulator.atmosphere.seeing['fwhm_ref'] = args.seeing * u.arcsec
+    simulator.atmosphere.seeing_fwhm_ref = args.seeing * u.arcsec
     seeing_fwhm = simulator.atmosphere.get_seeing_fwhm(
         wlen_grid).to(u.arcsec).value
 

specsim/transform.py

@@ -310,7 +310,10 @@ def create_observing_model(where, when, wavelength, temperature=15*u.deg_C,
     # See https://en.wikipedia.org/wiki/Vertical_pressure_variation
     if pressure is None:
         h = where.height
-        p0 = astropy.constants.atmosphere
+        if astropy.version.major == 1:
+            p0 = astropy.constants.atmosphere
+        else:
+            p0 = astropy.constants.atm
         g0 = astropy.constants.g0
         R = astropy.constants.R
         air_molar_mass = 0.0289644 * u.kg / u.mol