diff --git a/specsim/fiberloss.py b/specsim/fiberloss.py index 478d42b..77f1455 100644 --- a/specsim/fiberloss.py +++ b/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,13 +52,17 @@ 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( @@ -62,7 +70,6 @@ def __init__(self, fiber_diameter, wlen_grid, num_pixels=16, 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( diff --git a/specsim/quickfiberloss.py b/specsim/quickfiberloss.py index f7ccf4a..3621c16 100644 --- a/specsim/quickfiberloss.py +++ b/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 diff --git a/specsim/transform.py b/specsim/transform.py index fe926cf..83ff00b 100644 --- a/specsim/transform.py +++ b/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