From c6c6021d0750567d5589c29c37bb81937497fddb Mon Sep 17 00:00:00 2001 From: Luca Pagano Date: Tue, 5 Nov 2024 15:58:37 +0100 Subject: [PATCH] reverted the pointing API to the old behaviour, documentation updated --- docs/source/scanning.rst | 35 +++++++++++++++++++++++--------- litebird_sim/dipole.py | 2 +- litebird_sim/hwp_sys/hwp_sys.py | 1 - litebird_sim/mapmaking/common.py | 1 - litebird_sim/observations.py | 15 +++++++++++--- litebird_sim/scan_map.py | 1 - test/test_coordinates.py | 10 ++++----- test/test_pointing_sys.py | 12 +++++------ test/test_simulations.py | 12 +++++------ 9 files changed, 55 insertions(+), 34 deletions(-) diff --git a/docs/source/scanning.rst b/docs/source/scanning.rst index 5dc99e5d..91fa9eb5 100644 --- a/docs/source/scanning.rst +++ b/docs/source/scanning.rst @@ -168,15 +168,15 @@ similar to what is going to be used for LiteBIRD: .. testoutput:: - Shape: (1, 600, 3) + Shape: (600, 3) Pointings: - [[[ 2.182 -0. -1.571] + [[ 2.182 -0. -1.571] [ 2.182 -0.006 -1.576] [ 2.182 -0.012 -1.582] ... [ 0.089 -2.967 -1.738] [ 0.088 -3.021 -1.687] - [ 0.087 -3.075 -1.635]]] + [ 0.087 -3.075 -1.635]] All the details in this code are explained in the next sections, so for now just keep in mind the overall shape of the code: @@ -201,12 +201,27 @@ for now just keep in mind the overall shape of the code: the example above, this is done by the function :func:`.get_pointings`. -3. The method :meth:`.Observation.get_pointings` returns a ``(D, N, 3)`` - matrix, where D represents the detector index, N the index of the sample - and the three final columns contain the colatitude :math:`\theta`, - the longitude :math:`\phi`, and the orientation angle :math:`\psi`, - all expressed in radians. These angles are expressed in the Ecliptic - Coordinate System, where the Equator is aligned with the Ecliptic Plane of +3. The method :meth:`.Observation.get_pointings` returns an array with + either 2 or 3 fields depending on the argument passed: + + - if an integer is passed, this is interpreted as the index of the + detector in the observation, and a ``(N, 3)`` matrix is returned + where the first column contains the colatitude :math:`\theta`, + the second column the longitude :math:`\phi`, and the third column + the orientation angle :math:`\psi`, all expressed in radians. + + - if a list containing indices is passed, this is interpreted as + a list of detectors in the observation for which we want to compute + the pointing. It returns a ``(D, N, 3)`` matrix where D represents + the detector index, N the index of the sample and the three final + columns are the same described in the first case. + + - if the string "all" is passed then a ``(D, N, 3)`` matrix is returned + containig the pointing information for all the detectors in the + observation. + + These angles are expressed in the Ecliptic Coordinate + System, where the Equator is aligned with the Ecliptic Plane of the Solar System. @@ -361,7 +376,7 @@ split in several blocks inside the :class:`.Observation` class. unlikely to be relevant. Once all the quaternions have been computed at the proper sampling -rate, the direction of the detector on the sky and its o]rientation +rate, the direction of the detector on the sky and its orientation angle can be computed via a call to :meth:`.Observation.get_pointings`. The calculation works as follows: diff --git a/litebird_sim/dipole.py b/litebird_sim/dipole.py index ad0953ae..0ae381c0 100644 --- a/litebird_sim/dipole.py +++ b/litebird_sim/dipole.py @@ -225,7 +225,7 @@ def add_dipole( if type(pointings) is np.ndarray: theta_phi_det = pointings[detector_idx, :, :] else: - theta_phi_det = pointings(detector_idx)[0][0, :, 0:2] + theta_phi_det = pointings(detector_idx)[0][:, 0:2] add_dipole_for_one_detector( tod_det=tod[detector_idx], diff --git a/litebird_sim/hwp_sys/hwp_sys.py b/litebird_sim/hwp_sys/hwp_sys.py index d9b55f8d..0b6034e3 100644 --- a/litebird_sim/hwp_sys/hwp_sys.py +++ b/litebird_sim/hwp_sys/hwp_sys.py @@ -1256,7 +1256,6 @@ def fill_tod( cur_point, cur_hwp_angle = cur_obs.get_pointings( detector_idx=idet, pointings_dtype=dtype_pointings ) - cur_point = cur_point.reshape(-1, 3) else: cur_point = ptg_list[idx_obs][idet, :, :] cur_hwp_angle = hwp_angle_list[idx_obs] diff --git a/litebird_sim/mapmaking/common.py b/litebird_sim/mapmaking/common.py index fbd90699..d36929e9 100644 --- a/litebird_sim/mapmaking/common.py +++ b/litebird_sim/mapmaking/common.py @@ -205,7 +205,6 @@ def _compute_pixel_indices( curr_pointings_det = pointings[idet, :, :] else: curr_pointings_det, hwp_angle = pointings(idet) - curr_pointings_det = curr_pointings_det.reshape(-1, 3) if hwp_angle is None: hwp_angle = 0 diff --git a/litebird_sim/observations.py b/litebird_sim/observations.py index 7997de56..309e3101 100644 --- a/litebird_sim/observations.py +++ b/litebird_sim/observations.py @@ -680,12 +680,21 @@ def get_pointings( # All the detectors are included pointings, hwp_angle = cur_obs.get_pointings("all") + # This returns ``(N_det, N_samples, 3)`` array # Only the first two detectors are included pointings, hwp_angle = cur_obs.get_pointings([0, 1]) + # This returns ``(2, N_samples, 3)`` array # Only the first detector is used pointings, hwp_angle = cur_obs.get_pointings(0) + # This returns ``(N_samples, 3)`` array + + # NB if a list of indices is passed an array of dimension ``(N_det, N_samples, 3)`` + # is always returned + # For example this + pointings, hwp_angle = cur_obs.get_pointings([0]) + # returns ``(1, N_samples, 3)`` array The return value is a pair containing (1) the pointing matrix and (2) the HWP angle. The pointing matrix is a NumPy array with shape ``(N_det, N_samples, 3)``, @@ -693,7 +702,8 @@ def get_pointings( samples in the TOD (the field ``Observation.n_samples``). The last dimension spans the three angles θ (colatitude, in radians), φ (longitude, in radians), and ψ (orientation angle, in radians). *Important*: if you ask for just *one* - detector, the shape of the pointing matrix will always be ``(N_samples, 3)``. + detector passing the index of the detector, the shape of the pointing matrix + will always be ``(N_samples, 3)``. The HWP angle is always a vector with shape ``(N_samples,)``, as it does not depend on the list of detectors. @@ -712,7 +722,7 @@ def get_pointings( (detector_idx >= 0) and (detector_idx < self.n_detectors) ), f"Invalid detector index {detector_idx}, it must be a number between 0 and {self.n_detectors - 1}" - pointing, hwp = self.pointing_provider.get_pointings( + return self.pointing_provider.get_pointings( detector_quat=self.quat[detector_idx], start_time=self.start_time, start_time_global=self.start_time_global, @@ -722,7 +732,6 @@ def get_pointings( hwp_buffer=hwp_buffer, pointings_dtype=pointings_dtype, ) - return pointing.reshape(1, -1, 3), hwp # More complex case: we need all the detectors if isinstance(detector_idx, str): diff --git a/litebird_sim/scan_map.py b/litebird_sim/scan_map.py index 205c6ad1..23e24c3d 100644 --- a/litebird_sim/scan_map.py +++ b/litebird_sim/scan_map.py @@ -63,7 +63,6 @@ def scan_map( curr_pointings_det = pointings[detector_idx, :, :] else: curr_pointings_det, hwp_angle = pointings(detector_idx) - curr_pointings_det = curr_pointings_det.reshape(-1, 3) if hwp_angle is None: hwp_angle = 0 diff --git a/test/test_coordinates.py b/test/test_coordinates.py index 2ab08288..b7ebb983 100644 --- a/test/test_coordinates.py +++ b/test/test_coordinates.py @@ -49,14 +49,14 @@ def test_coordinates(): r = hp.Rotator(coord=["E", "G"]) - pointings_gal_hp = np.empty_like(pointings[0]) + pointings_gal_hp = np.empty_like(pointings) - pointings_gal_hp[:, 0:2] = r(pointings[0, :, 0], pointings[0, :, 1]).T - pointings_gal_hp[:, 2] = pointings[0, :, 2] + r.angle_ref( - pointings[0, :, 0], pointings[0, :, 1] + pointings_gal_hp[:, 0:2] = r(pointings[:, 0], pointings[:, 1]).T + pointings_gal_hp[:, 2] = pointings[:, 2] + r.angle_ref( + pointings[:, 0], pointings[:, 1] ) - pointings_gal_lbs = lbs.coordinates.rotate_coordinates_e2g(pointings[0]) + pointings_gal_lbs = lbs.coordinates.rotate_coordinates_e2g(pointings) np.testing.assert_allclose( pointings_gal_hp, pointings_gal_lbs, rtol=1e-6, atol=1e-6 diff --git a/test/test_pointing_sys.py b/test/test_pointing_sys.py index faaeafe3..c89cebf2 100644 --- a/test/test_pointing_sys.py +++ b/test/test_pointing_sys.py @@ -112,7 +112,7 @@ def test_PointingSys_add_single_offset_to_FP( pointings, hwp_angle = cur_obs.get_pointings( det_idx, pointings_dtype=np.float32 ) - pointings_list.append(pointings[0].tolist()) + pointings_list.append(pointings.tolist()) if func.__name__ not in results_dict: results_dict[func.__name__] = {} @@ -150,7 +150,7 @@ def test_PointingSys_add_multiple_offsets_to_FP( pointings, hwp_angle = cur_obs.get_pointings( det_idx, pointings_dtype=np.float32 ) - pointings_list.append(pointings[0].tolist()) + pointings_list.append(pointings.tolist()) if func.__name__ not in results_dict: results_dict[func.__name__] = {} @@ -196,7 +196,7 @@ def test_PointingSys_add_uncommon_disturb_to_FP( pointings, hwp_angle = cur_obs.get_pointings( det_idx, pointings_dtype=np.float32 ) - pointings_list.append(pointings[0].tolist()) + pointings_list.append(pointings.tolist()) if func.__name__ not in results_dict: results_dict[func.__name__] = {} @@ -241,7 +241,7 @@ def test_PointingSys_add_common_disturb_to_FP( pointings, hwp_angle = cur_obs.get_pointings( det_idx, pointings_dtype=np.float32 ) - pointings_list.append(pointings[0].tolist()) + pointings_list.append(pointings.tolist()) if func.__name__ not in results_dict: results_dict[func.__name__] = {} @@ -279,7 +279,7 @@ def test_PointingSys_add_single_offset_to_spacecraft( pointings, hwp_angle = cur_obs.get_pointings( det_idx, pointings_dtype=np.float32 ) - pointings_list.append(pointings[0].tolist()) + pointings_list.append(pointings.tolist()) if func.__name__ not in results_dict: results_dict[func.__name__] = {} @@ -323,7 +323,7 @@ def test_PointingSys_add_common_disturb_to_spacecraft( pointings, hwp_angle = cur_obs.get_pointings( det_idx, pointings_dtype=np.float32 ) - pointings_list.append(pointings[0].tolist()) + pointings_list.append(pointings.tolist()) if func.__name__ not in results_dict: results_dict[func.__name__] = {} diff --git a/test/test_simulations.py b/test/test_simulations.py index 4b77f020..fb5b8a57 100644 --- a/test/test_simulations.py +++ b/test/test_simulations.py @@ -465,7 +465,7 @@ def test_smart_pointings_consistency_with_hwp(tmp_path): for det_idx in range(obs.n_detectors): (pointings, hwp_angle) = obs.get_pointings(det_idx) - assert pointings.shape == (1, obs.n_samples, 3) + assert pointings.shape == (obs.n_samples, 3) assert hwp_angle.shape == (obs.n_samples,) @@ -479,7 +479,7 @@ def test_smart_pointings_consistency_without_hwp(tmp_path): for det_idx in range(obs.n_detectors): (pointings, hwp_angle) = obs.get_pointings(det_idx) - assert pointings.shape == (1, obs.n_samples, 3) + assert pointings.shape == (obs.n_samples, 3) assert hwp_angle is None @@ -497,7 +497,7 @@ def test_smart_pointings_angles(tmp_path): # function present in LBS 0.12.0 np.testing.assert_allclose( - actual=pointings[0, 0:10, :], + actual=pointings[0:10, :], desired=np.array( [ [1.6580627894, 0.0000000000, 1.3977241805], @@ -617,7 +617,7 @@ def test_compute_pointings_for_one_detector(dtype, tmp_path): pointings, hwp_angle = cur_obs.get_pointings(0, pointings_dtype=dtype) assert pointings.dtype == dtype - assert pointings.shape == (1, cur_obs.n_samples, 3) + assert pointings.shape == (cur_obs.n_samples, 3) assert hwp_angle.dtype == dtype assert hwp_angle.shape == (cur_obs.n_samples,) @@ -648,7 +648,7 @@ def test_store_pointings_for_two_detectors(dtype, tmp_path): cur_pointings, _ = cur_obs.get_pointings(abs_det_idx) np.testing.assert_allclose( actual=pointings[rel_det_idx, :, :], - desired=cur_pointings[0, :, :], + desired=cur_pointings[:, :], rtol=1e-6, ) @@ -677,6 +677,6 @@ def test_smart_pointings_for_all_detectors(dtype, tmp_path): assert cur_pointings.dtype == dtype np.testing.assert_allclose( actual=pointings[det_idx, :, :], - desired=cur_pointings[0, :, :], + desired=cur_pointings[:, :], rtol=1e-6, )