diff --git a/scopesim/effects/electronic/noise.py b/scopesim/effects/electronic/noise.py index 3946b4dd..cb363412 100644 --- a/scopesim/effects/electronic/noise.py +++ b/scopesim/effects/electronic/noise.py @@ -136,32 +136,33 @@ def __init__(self, **kwargs): self.meta.update(kwargs) def apply_to(self, det, **kwargs): - if isinstance(det, DetectorBase): - self.meta["random_seed"] = from_currsys(self.meta["random_seed"], - self.cmds) - if self.meta["random_seed"] is not None: - np.random.seed(self.meta["random_seed"]) - - # ! poisson(x) === normal(mu=x, sigma=x**0.5) - # Windows has a problem with generating poisson values above 2**30 - # Above ~100 counts the poisson and normal distribution are - # basically the same. For large arrays the normal distribution - # takes only 60% as long as the poisson distribution - data = det._hdu.data - - # Check if there are negative values in the data - negvals_mask = data < 0 - if negvals_mask.any(): - logger.warning(f"Effect ShotNoise: {negvals_mask.sum()} negative pixels") - data[negvals_mask] = 0 - - below = data < 2**20 - above = np.invert(below) - data[below] = np.random.poisson(data[below]).astype(float) - data[above] = np.random.normal(data[above], np.sqrt(data[above])) - new_imagehdu = fits.ImageHDU(data=data, header=det._hdu.header) - det._hdu = new_imagehdu - + if not isinstance(det, DetectorBase): + return det + + self.meta["random_seed"] = from_currsys(self.meta["random_seed"], + self.cmds) + rng = np.random.default_rng(self.meta["random_seed"]) + + # ! poisson(x) === normal(mu=x, sigma=x**0.5) + # Windows has a problem with generating poisson values above 2**30 + # Above ~100 counts the poisson and normal distribution are + # basically the same. For large arrays the normal distribution + # takes only 60% as long as the poisson distribution + + # Check if there are negative values in the data + data = np.ma.masked_less(det._hdu.data, 0) + if data.mask.any(): + logger.warning( + "Effect ShotNoise: %d negative pixels", data.mask.sum()) + data = data.filled(0) + + # Weirdly, poisson doesn't understand masked arrays, but normal does... + data = np.ma.masked_greater(data, 2**20) + data[~data.mask] = rng.poisson(data[~data.mask].data) + data.mask = ~data.mask + new_imagehdu = fits.ImageHDU(data=rng.normal(data, np.sqrt(data)), + header=det._hdu.header) + det._hdu = new_imagehdu return det def plot(self, det): diff --git a/scopesim/effects/obs_strategies.py b/scopesim/effects/obs_strategies.py index b34bd37e..a7b76aad 100644 --- a/scopesim/effects/obs_strategies.py +++ b/scopesim/effects/obs_strategies.py @@ -1,3 +1,4 @@ +# -*- coding: utf-8 -*- """ Effects describing observing strategies. @@ -67,21 +68,21 @@ def __init__(self, **kwargs): self.meta.update(kwargs) def apply_to(self, obj, **kwargs): - if isinstance(obj, DetectorBase): - chop_offsets = from_currsys(self.meta["chop_offsets"], self.cmds) - nod_offsets = from_currsys(self.meta["nod_offsets"], self.cmds) - if nod_offsets is None: - nod_offsets = -np.array(chop_offsets) - - # these offsets are in pixels, not in arcsec or mm - pixel_scale = float(from_currsys(self.meta["pixel_scale"], self.cmds)) - chop_offsets_pixel = np.array(chop_offsets) / pixel_scale - nod_offsets_pixel = np.array(nod_offsets) / pixel_scale - - image = obj.hdu.data - obj.hdu.data = chop_nod_image(image, - chop_offsets_pixel.astype(int), - nod_offsets_pixel.astype(int)) + if not isinstance(obj, DetectorBase): + return obj + + chop_offsets = from_currsys(self.meta["chop_offsets"], self.cmds) + nod_offsets = from_currsys(self.meta["nod_offsets"], self.cmds) + if nod_offsets is None: + nod_offsets = -np.array(chop_offsets) + + # these offsets are in pixels, not in arcsec or mm + pixel_scale = float(from_currsys(self.meta["pixel_scale"], self.cmds)) + chop_offsets_px = (np.array(chop_offsets) / pixel_scale).astype(int) + nod_offsets_px = (np.array(nod_offsets) / pixel_scale).astype(int) + + image = obj.hdu.data + obj.hdu.data = chop_nod_image(image, chop_offsets_px, nod_offsets_px) return obj @@ -97,5 +98,4 @@ def chop_nod_image(img, chop_offsets, nod_offsets=None): im_bb = np.roll(im_ba, chop_offsets, (1, 0)) im_comb = (im_aa - im_ab) - (im_ba - im_bb) - return im_comb diff --git a/scopesim/tests/tests_effects/test_ReferencePixelBorder.py b/scopesim/tests/tests_effects/test_ReferencePixelBorder.py index 7196dba1..e9edcf23 100644 --- a/scopesim/tests/tests_effects/test_ReferencePixelBorder.py +++ b/scopesim/tests/tests_effects/test_ReferencePixelBorder.py @@ -1,13 +1,8 @@ -import pytest -from pytest import approx - import numpy as np from matplotlib import pyplot as plt -from astropy import units as u from scopesim.effects import electronic as ee from scopesim.optics.image_plane import ImagePlane -from scopesim.base_classes import ImagePlaneBase from scopesim.tests.mocks.py_objects.imagehdu_objects import _image_hdu_square @@ -18,18 +13,18 @@ class TestInit: def test_initialised_with_nothing(self): rpb = ee.ReferencePixelBorder() assert isinstance(rpb, ee.ReferencePixelBorder) - assert np.all([rpb.meta[key] == 0 - for key in ["top", "bottom", "right", "left"]]) + assert all(rpb.meta[key] == 0 + for key in ["top", "bottom", "right", "left"]) def test_borders_all_set_to_5_for_keyword_all(self): rpb = ee.ReferencePixelBorder(all=5) - assert np.all([rpb.meta[key] == 5 - for key in ["top", "bottom", "right", "left"]]) + assert all(rpb.meta[key] == 5 + for key in ["top", "bottom", "right", "left"]) def test_border_set_differently(self): rpb = ee.ReferencePixelBorder(top=5, right=3) borders = {"top": 5, "bottom": 0, "right": 3, "left": 0} - assert np.all([rpb.meta[key] == borders[key] for key in borders]) + assert all(rpb.meta[key] == borders[key] for key in borders) class TestApplyTo: diff --git a/scopesim/tests/tests_effects/test_fits_headers.py b/scopesim/tests/tests_effects/test_fits_headers.py index 21c3197a..912ba200 100644 --- a/scopesim/tests/tests_effects/test_fits_headers.py +++ b/scopesim/tests/tests_effects/test_fits_headers.py @@ -145,7 +145,7 @@ def test_works_for_ext_name(self, comb_hdul): exts = fh.get_relevant_extensions(dic, comb_hdul) answer = [0] - assert np.all([ans in exts for ans in answer]) + assert all(ans in exts for ans in answer) assert len(exts) == len(answer) def test_works_for_ext_number(self, comb_hdul): @@ -153,7 +153,7 @@ def test_works_for_ext_number(self, comb_hdul): exts = fh.get_relevant_extensions(dic, comb_hdul) answer = [1, 2] - assert np.all([ans in exts for ans in answer]) + assert all(ans in exts for ans in answer) assert len(exts) == len(answer) @pytest.mark.parametrize("ext_type, answer", @@ -163,7 +163,7 @@ def test_works_for_ext_type(self, comb_hdul, ext_type, answer): dic = {"ext_type": ext_type} exts = fh.get_relevant_extensions(dic, comb_hdul) - assert np.all([ans in exts for ans in answer]) + assert all(ans in exts for ans in answer) assert len(exts) == len(answer) diff --git a/scopesim/tests/tests_effects/test_obs_strategies.py b/scopesim/tests/tests_effects/test_obs_strategies.py index bcc4489e..47001c54 100644 --- a/scopesim/tests/tests_effects/test_obs_strategies.py +++ b/scopesim/tests/tests_effects/test_obs_strategies.py @@ -1,9 +1,7 @@ import pytest -from pytest import approx import numpy as np import matplotlib.pyplot as plt -from matplotlib.colors import LogNorm from scopesim.effects.obs_strategies import ChopNodCombiner from scopesim.detector import Detector @@ -44,7 +42,9 @@ def test_creates_image_for_parallel_chop(self, basic_detector): plt.imshow(basic_detector.hdu.data) plt.show() - assert np.sum(basic_detector.hdu.data) == 0 + outimg = basic_detector.hdu.data + assert outimg.sum() == 0 + assert ((outimg == 0.).sum() / outimg.size) > .8 # most elements zero def test_creates_image_for_perpendicular_chop(self, basic_detector): cnc = ChopNodCombiner(pixel_scale=0.004, chop_offsets=(0.12, 0), @@ -55,4 +55,6 @@ def test_creates_image_for_perpendicular_chop(self, basic_detector): plt.imshow(basic_detector.hdu.data) plt.show() - assert np.sum(basic_detector.hdu.data) == 0 + outimg = basic_detector.hdu.data + assert outimg.sum() == 0 + assert ((outimg == 0.).sum() / outimg.size) > .8 # most elements zero diff --git a/scopesim/tests/tests_optics/test_FieldOfView.py b/scopesim/tests/tests_optics/test_FieldOfView.py index b4b26d35..bf27ef04 100644 --- a/scopesim/tests/tests_optics/test_FieldOfView.py +++ b/scopesim/tests/tests_optics/test_FieldOfView.py @@ -3,6 +3,7 @@ import numpy as np from astropy import units as u from astropy.io import fits +from astropy.table import Table from matplotlib import pyplot as plt from matplotlib.colors import LogNorm @@ -129,21 +130,19 @@ def test_all_spectra_are_referenced_correctly(self): # check the same spectrum object is referenced by both lists assert fov.fields[0].header["SPEC_REF"] == \ src.fields[0].header["SPEC_REF"] - assert np.all([fov.fields[2][i]["ref"] == \ - src.fields[2][i]["ref"] for i in range(4)]) - - def test_contains_all_fields_inside_fov(self, basic_fov_header, - cube_source, - image_source, table_source): - src = image_source + cube_source + table_source - the_fov = FieldOfView(basic_fov_header, (1, 2) * u.um, - area=1 * u.m ** 2) - the_fov.extract_from(src) - assert len(the_fov.fields) == 3 - assert isinstance(the_fov.fields[0], fits.ImageHDU) - assert isinstance(the_fov.fields[1], fits.ImageHDU) - assert the_fov.fields[1].header["NAXIS"] == 3 - assert isinstance(the_fov.fields[2], Table) + assert all(fov.fields[2][i]["ref"] == src.fields[2][i]["ref"] + for i in range(4)) + + def test_contains_all_fields_inside_fov(self): + src = so._image_source() + so._cube_source() + so._table_source() + the_fov = FieldOfView(ho._basic_fov_header(), (1, 2) * u.um, + area=1 * u.m ** 2) + the_fov.extract_from(src) + assert len(the_fov.fields) == 3 + assert isinstance(the_fov.fields[0], fits.ImageHDU) + assert isinstance(the_fov.fields[1], fits.ImageHDU) + assert the_fov.fields[1].header["NAXIS"] == 3 + assert isinstance(the_fov.fields[2], Table) def test_handles_nans(self): src = so._image_source() diff --git a/scopesim/tests/tests_optics/test_SpectralSurface.py b/scopesim/tests/tests_optics/test_SpectralSurface.py index 17c83485..ac17234b 100644 --- a/scopesim/tests/tests_optics/test_SpectralSurface.py +++ b/scopesim/tests/tests_optics/test_SpectralSurface.py @@ -153,7 +153,7 @@ def test_the_right_answers_for_valid_input(self, colname1, colname2, srf.meta[colname1] = col1 srf.meta[colname2] = col2 col3 = srf._compliment_array(colname1, colname2) - assert np.all(np.isclose(col3.data, expected.data)) + assert np.allclose(col3.data, expected.data) assert col3.unit == expected.unit @pytest.mark.parametrize("colname1, colname2, col1, col2, expected", diff --git a/scopesim/tests/tests_server/test_database.py b/scopesim/tests/tests_server/test_database.py index 07b15900..9491cb06 100644 --- a/scopesim/tests/tests_server/test_database.py +++ b/scopesim/tests/tests_server/test_database.py @@ -4,7 +4,6 @@ from tempfile import TemporaryDirectory import yaml -import numpy as np from scopesim.server import database as db from scopesim.server import example_data_utils as dbex @@ -95,7 +94,7 @@ def test_lists_all_packages_without_qualifier(self): def test_lists_only_packages_with_qualifier(self): pkgs = db.list_packages("Armazones") - assert np.all(["Armazones" in pkg for pkg in pkgs]) + assert all("Armazones" in pkg for pkg in pkgs) def test_throws_for_nonexisting_pkgname(self): with pytest.raises(ValueError): diff --git a/scopesim/tests/tests_source/test_source_Source.py b/scopesim/tests/tests_source/test_source_Source.py index 3859ee5d..0d8dda20 100644 --- a/scopesim/tests/tests_source/test_source_Source.py +++ b/scopesim/tests/tests_source/test_source_Source.py @@ -327,14 +327,16 @@ def test_photons_returned_only_for_indices(self, table_source): assert np.allclose(ph.value, [4, 2]) +@pytest.mark.xfail class TestSourceShift: def test_that_it_does_what_it_should(self): - pass + assert False +@pytest.mark.xfail class TestSourceRotate: def test_that_it_does_what_it_should(self): - pass + assert False class TestPhotonsInRange: