From 909ce8c19dc3aeb5ffe6fdc864490cb34f0190d8 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Wed, 18 Apr 2018 12:26:44 +0200 Subject: [PATCH 01/14] Solving crop_by_coord problem when data are rotated --- ndcube/ndcube.py | 75 +++++++++++++++++++++++++++++++++++------------- 1 file changed, 55 insertions(+), 20 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index 696e7b3b0..bcf5fcb9e 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -1,6 +1,8 @@ # -*- coding: utf-8 -*- import abc +import warnings +from itertools import product import numpy as np import astropy.nddata @@ -85,22 +87,30 @@ def world_axis_physical_types(self): pass @abc.abstractmethod - def crop_by_coords(self, min_coord_values, interval_widths): + def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): """ Crops an NDCube given minimum values and interval widths along axes. Parameters ---------- - min_coord_values: iterable of `astropy.units.Quantity` + lower_corner: iterable of `astropy.units.Quantity` The minimum desired values along each relevant axis after cropping described in physical units consistent with the NDCube's wcs object. - The length of the iterable must equal the number of data dimensions and must - have the same order as the data. + The length of the iterable must equal the number of data dimensions + and must have the same order as the data. interval_widths: iterable of `astropy.units.Quantity` - The width of the region of interest in each dimension in physical units - consistent with the NDCube's wcs object. The length of the iterable must - equal the number of data dimensions and must have the same order as the data. + The width of the region of interest in each dimension in physical + units consistent with the NDCube's wcs object. The length of the + iterable must equal the number of data dimensions and must have + the same order as the data. This argument will be removed in versions + 2.0, please use upper_corner argument. + + upper_corner: iterable of `astropy.units.Quantity` + The maximum desired values along each relevant axis after cropping + described in physical units consistent with the NDCube's wcs object. + The length of the iterable must equal the number of data dimensions + and must have the same order as the data. Returns ------- @@ -409,21 +419,46 @@ def extra_coords(self): "value": self._extra_coords_wcs_axis[key]["value"]} return result - def crop_by_coords(self, min_coord_values, interval_widths): + def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): # The docstring is defined in NDDataBase n_dim = len(self.dimensions) - if (len(min_coord_values) != len(interval_widths)) or len(min_coord_values) != n_dim: - raise ValueError("min_coord_values and interval_widths must have " - "same number of elements as number of data dimensions.") - # Convert coords of lower left corner to pixel units. - lower_pixels = self.world_to_pixel(*min_coord_values) - upper_pixels = self.world_to_pixel(*[min_coord_values[i] + interval_widths[i] - for i in range(n_dim)]) - # Round pixel values to nearest integer. - lower_pixels = [int(np.rint(l.value)) for l in lower_pixels] - upper_pixels = [int(np.rint(u.value)) for u in upper_pixels] - item = tuple([slice(lower_pixels[i], upper_pixels[i]) for i in range(n_dim)]) + # Raising a value error if the arguments have not the same dimensions. + if upper_corner: + if (len(lower_corner) != len(upper_corner)) or (len(lower_corner) != n_dim): + raise ValueError("lower_corner and upper_corner must have " + "same number of elements as number of data " + "dimensions.") + # Raising a value error if the arguments have not the same dimensions. + # Calculation of upper_corner with the inputing interval_widths + # This part of the code will be removed in version 2.0 + if interval_widths: + warnings.warn("interval_widths will be removed from the API in " + "version 2.0, please use upper_corner argument.") + if (len(lower_corner) != len(interval_widths)) or (len(lower_corner) != n_dim): + raise ValueError("lower_corner and interval_widths must have " + "same number of elements as number of data " + "dimensions.") + upper_corner = [lower_corner[i] + interval_widths[i] for i in range(n_dim)] + # Derive all corners coordinates + quantity_list = [[lower_corner[i], upper_corner[i]] for i in range(n_dim)] + all_corners = [self.world_to_pixel(*a) for a in product(*quantity_list)] + # Inputing of all corners coordinates inside a numpy array + # According to the boundary conditions + corners_array = np.zeros((2**n_dim, n_dim)) + for i in range(2**n_dim): + for j in range(n_dim): + corners_array[i, j] = u.Quantity(all_corners[i][j]).value + if corners_array[i, j] > self.data.shape[j]: + corners_array[i, j] = self.data.shape[j] + if corners_array[i, j] < 0: + corners_array[i, j] = 0 + # Taking the maximum and minimum values of coordinates + lower_pixels = corners_array.min(0) + upper_pixels = corners_array.max(0) + # Creating a tuple to crop the data with inputed coordinates + item = tuple([slice(int(lower_pixels[i]), int(upper_pixels[i])) for i in range(n_dim)]) + return self[item] def crop_by_extra_coord(self, min_coord_value, interval_width, coord_name): @@ -502,7 +537,7 @@ class NDCubeOrdered(NDCube): for standard deviation or "var" for variance. A metaclass defining such an interface is NDUncertainty - but isn’t mandatory. If the uncertainty has no such attribute the uncertainty is stored as UnknownUncertainty. - Defaults to None. + Defaults to None.list mask : any type, optional Mask for the dataset. Masks should follow the numpy convention From 28aba8cfa5bc02ae19d630b2a8ebd132dae8dffc Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Thu, 19 Apr 2018 13:57:12 +0200 Subject: [PATCH 02/14] New error raising and dimension modified --- ndcube/ndcube.py | 25 ++++++++++++++----------- 1 file changed, 14 insertions(+), 11 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index bcf5fcb9e..08c2add90 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -422,24 +422,27 @@ def extra_coords(self): def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): # The docstring is defined in NDDataBase - n_dim = len(self.dimensions) - # Raising a value error if the arguments have not the same dimensions. - if upper_corner: - if (len(lower_corner) != len(upper_corner)) or (len(lower_corner) != n_dim): - raise ValueError("lower_corner and upper_corner must have " - "same number of elements as number of data " - "dimensions.") + n_dim = self.data.ndim # Raising a value error if the arguments have not the same dimensions. # Calculation of upper_corner with the inputing interval_widths # This part of the code will be removed in version 2.0 if interval_widths: warnings.warn("interval_widths will be removed from the API in " "version 2.0, please use upper_corner argument.") + if upper_corner: + raise ValueError("Only one of interval_widths or upper_corner " + "can be set. Recommend using upper_corner as " + "interval_widths is deprecated.") if (len(lower_corner) != len(interval_widths)) or (len(lower_corner) != n_dim): raise ValueError("lower_corner and interval_widths must have " "same number of elements as number of data " "dimensions.") - upper_corner = [lower_corner[i] + interval_widths[i] for i in range(n_dim)] + upper_corner = [lower_corner[i] + interval_widths[i] + for i in range(n_dim)] + # Raising a value error if the arguments have not the same dimensions. + if (len(lower_corner) != len(upper_corner)) or (len(lower_corner) != n_dim): + raise ValueError("lower_corner and upper_corner must have same" + "number of elements as number of data dimensions.") # Derive all corners coordinates quantity_list = [[lower_corner[i], upper_corner[i]] for i in range(n_dim)] all_corners = [self.world_to_pixel(*a) for a in product(*quantity_list)] @@ -457,8 +460,8 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): lower_pixels = corners_array.min(0) upper_pixels = corners_array.max(0) # Creating a tuple to crop the data with inputed coordinates - item = tuple([slice(int(lower_pixels[i]), int(upper_pixels[i])) for i in range(n_dim)]) - + item = tuple([slice(int(lower_pixels[i]), int(upper_pixels[i])) + for i in range(n_dim)]) return self[item] def crop_by_extra_coord(self, min_coord_value, interval_width, coord_name): @@ -537,7 +540,7 @@ class NDCubeOrdered(NDCube): for standard deviation or "var" for variance. A metaclass defining such an interface is NDUncertainty - but isn’t mandatory. If the uncertainty has no such attribute the uncertainty is stored as UnknownUncertainty. - Defaults to None.list + Defaults to None. mask : any type, optional Mask for the dataset. Masks should follow the numpy convention From 952a71296742bafaaa92324434601c656ade03b8 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Fri, 20 Apr 2018 15:37:17 +0200 Subject: [PATCH 03/14] Adding a new test for rotated data --- ndcube/ndcube.py | 6 ++---- ndcube/tests/test_ndcube.py | 27 ++++++++++++++++++++++++++- 2 files changed, 28 insertions(+), 5 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index 08c2add90..d169b5e2d 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -437,8 +437,7 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): raise ValueError("lower_corner and interval_widths must have " "same number of elements as number of data " "dimensions.") - upper_corner = [lower_corner[i] + interval_widths[i] - for i in range(n_dim)] + upper_corner = [lower_corner[i] + interval_widths[i] for i in range(n_dim)] # Raising a value error if the arguments have not the same dimensions. if (len(lower_corner) != len(upper_corner)) or (len(lower_corner) != n_dim): raise ValueError("lower_corner and upper_corner must have same" @@ -460,8 +459,7 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): lower_pixels = corners_array.min(0) upper_pixels = corners_array.max(0) # Creating a tuple to crop the data with inputed coordinates - item = tuple([slice(int(lower_pixels[i]), int(upper_pixels[i])) - for i in range(n_dim)]) + item = tuple([slice(int(lower_pixels[i]), int(upper_pixels[i]) + 1) for i in range(n_dim)]) return self[item] def crop_by_extra_coord(self, min_coord_value, interval_width, coord_name): diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index b400d1531..151c3220a 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -59,6 +59,20 @@ data_ordered[0] = data.transpose() data_ordered[1] = data.transpose() +h_rotated = {'CTYPE1' : 'HPLN-TAN', 'CUNIT1' : 'arcsec', 'CDELT1' : 0.4, 'CRPIX1' : 0, + 'CRVAL1' : 0, 'NAXIS1' : 5, + 'CTYPE2' : 'HPLT-TAN', 'CUNIT2' : 'arcsec', 'CDELT2' : 0.5, 'CRPIX2' : 0, + 'CRVAL2' : 0, 'NAXIS2' : 5, + 'CTYPE3' : 'Time ', 'CUNIT3' : 'seconds', 'CDELT3' : 0.3, 'CRPIX3' : 0, + 'CRVAL3' : 0, 'NAXIS3' : 2, + 'PC1_1' : 0.714963912964, 'PC1_2' : -0.699137151241, 'PC1_3' : 0.0, + 'PC2_1' : 0.699137151241, 'PC2_2' : 0.714963912964, 'PC2_3' : 0.0, + 'PC3_1' : 0.0, 'PC3_2' : 0.0, 'PC3_3' : 1.0} +w_rotated = WCS(header=h_rotated, naxis=3) + +data_rotated = np.array([[[1, 2, 3, 4, 6], [2, 4, 5, 3, 1], [0, -1, 2, 4, 2], [3, 5, 1, 2, 0]], + [[2, 4, 5, 1, 3], [1, 5, 2, 2, 4], [2, 3, 4, 0, 5], [0, 1, 2, 3, 4]]]) + mask_cubem = data > 0 mask_cube = data >= 0 uncertaintym = data @@ -117,6 +131,15 @@ ('hello', 1, u.Quantity(range(data.shape[1]), unit=u.pix)), ('bye', 2, u.Quantity(range(data.shape[2]), unit=u.pix))]) +cube_rotated = NDCube( + data_rotated, + w_rotated, + mask=mask_cube, + uncertainty=uncertainty, + missing_axis=[False, False, False], + extra_coords=[('time', 0, u.Quantity(range(data_rotated.shape[0]), unit=u.pix)), + ('hello', 1, u.Quantity(range(data_rotated.shape[1]), unit=u.pix)), + ('bye', 2,u.Quantity(range(data_rotated.shape[2]), unit=u.pix))]) @pytest.mark.parametrize( "test_input,expected,mask,wcs,uncertainty,dimensions,world_axis_physical_types,extra_coords", @@ -814,7 +837,9 @@ def test_world_to_pixel(test_input, expected): @pytest.mark.parametrize("test_input,expected", [ - ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 1.06*u.m]), cubem[:, :2])]) + ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 1.06*u.m]), cubem[:, :2]), + ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], [1*u.s, 1*u.arcsec, 0.5*u.arcsec]), + cube_rotated[:, :3, 1:4])]) def test_crop_by_coords(test_input, expected): helpers.assert_cubes_equal( test_input[0].crop_by_coords(*test_input[1:]), expected) From cefd9a7bfb08e751c95286fdf62daa1464989a28 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Fri, 20 Apr 2018 15:41:11 +0200 Subject: [PATCH 04/14] PEP8 correction in test_ndcube.py --- ndcube/tests/test_ndcube.py | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index 151c3220a..546c434cf 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -59,15 +59,15 @@ data_ordered[0] = data.transpose() data_ordered[1] = data.transpose() -h_rotated = {'CTYPE1' : 'HPLN-TAN', 'CUNIT1' : 'arcsec', 'CDELT1' : 0.4, 'CRPIX1' : 0, - 'CRVAL1' : 0, 'NAXIS1' : 5, - 'CTYPE2' : 'HPLT-TAN', 'CUNIT2' : 'arcsec', 'CDELT2' : 0.5, 'CRPIX2' : 0, - 'CRVAL2' : 0, 'NAXIS2' : 5, - 'CTYPE3' : 'Time ', 'CUNIT3' : 'seconds', 'CDELT3' : 0.3, 'CRPIX3' : 0, - 'CRVAL3' : 0, 'NAXIS3' : 2, - 'PC1_1' : 0.714963912964, 'PC1_2' : -0.699137151241, 'PC1_3' : 0.0, - 'PC2_1' : 0.699137151241, 'PC2_2' : 0.714963912964, 'PC2_3' : 0.0, - 'PC3_1' : 0.0, 'PC3_2' : 0.0, 'PC3_3' : 1.0} +h_rotated = {'CTYPE1': 'HPLN-TAN', 'CUNIT1': 'arcsec', 'CDELT1': 0.4, 'CRPIX1': 0, + 'CRVAL1': 0, 'NAXIS1': 5, + 'CTYPE2': 'HPLT-TAN', 'CUNIT2': 'arcsec', 'CDELT2': 0.5, 'CRPIX2': 0, + 'CRVAL2': 0, 'NAXIS2': 5, + 'CTYPE3': 'Time ', 'CUNIT3': 'seconds', 'CDELT3': 0.3, 'CRPIX3': 0, + 'CRVAL3': 0, 'NAXIS3': 2, + 'PC1_1': 0.714963912964, 'PC1_2': -0.699137151241, 'PC1_3': 0.0, + 'PC2_1': 0.699137151241, 'PC2_2': 0.714963912964, 'PC2_3': 0.0, + 'PC3_1': 0.0, 'PC3_2': 0.0, 'PC3_3': 1.0} w_rotated = WCS(header=h_rotated, naxis=3) data_rotated = np.array([[[1, 2, 3, 4, 6], [2, 4, 5, 3, 1], [0, -1, 2, 4, 2], [3, 5, 1, 2, 0]], @@ -139,7 +139,7 @@ missing_axis=[False, False, False], extra_coords=[('time', 0, u.Quantity(range(data_rotated.shape[0]), unit=u.pix)), ('hello', 1, u.Quantity(range(data_rotated.shape[1]), unit=u.pix)), - ('bye', 2,u.Quantity(range(data_rotated.shape[2]), unit=u.pix))]) + ('bye', 2, u.Quantity(range(data_rotated.shape[2]), unit=u.pix))]) @pytest.mark.parametrize( "test_input,expected,mask,wcs,uncertainty,dimensions,world_axis_physical_types,extra_coords", From 5a0f4540b0afc8e7075ddceb0837fb0da9ffd896 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Mon, 23 Apr 2018 11:32:33 +0200 Subject: [PATCH 05/14] Solving test bug with old values on crop by coords --- ndcube/tests/test_ndcube.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index 546c434cf..8c936827b 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -837,7 +837,7 @@ def test_world_to_pixel(test_input, expected): @pytest.mark.parametrize("test_input,expected", [ - ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 1.06*u.m]), cubem[:, :2]), + ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 5.e-11*u.m]), cubem[:, :, :3]), ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], [1*u.s, 1*u.arcsec, 0.5*u.arcsec]), cube_rotated[:, :3, 1:4])]) def test_crop_by_coords(test_input, expected): From b837d0ec5a4511fc113896403b240d97ff3bdb5f Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Mon, 23 Apr 2018 16:47:37 +0200 Subject: [PATCH 06/14] Some error tests added --- ndcube/tests/test_ndcube.py | 17 ++++++++++++----- 1 file changed, 12 insertions(+), 5 deletions(-) diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index 8c936827b..10b004dc2 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -837,8 +837,13 @@ def test_world_to_pixel(test_input, expected): @pytest.mark.parametrize("test_input,expected", [ - ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 5.e-11*u.m]), cubem[:, :, :3]), - ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], [1*u.s, 1*u.arcsec, 0.5*u.arcsec]), + ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 5.e-11*u.m], None), + cubem[:, :, :3]), + ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], [1*u.s, 1*u.arcsec, 0.5*u.arcsec], None), + cube_rotated[:, :3, 1:4]), + ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], None, [1.7*u.deg, 1*u.deg, 1.07e-9*u.m]), + cubem[:, :, :3]), + ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], None, [1*u.s, 2.5*u.arcsec, 0.5*u.arcsec]), cube_rotated[:, :3, 1:4])]) def test_crop_by_coords(test_input, expected): helpers.assert_cubes_equal( @@ -846,10 +851,12 @@ def test_crop_by_coords(test_input, expected): @pytest.mark.parametrize("test_input", [ - (cubem, u.Quantity([0], unit=u.deg), u.Quantity([1.5, 2.], unit=u.deg))]) + (ValueError, cubem, u.Quantity([0], unit=u.deg), u.Quantity([1.5, 2.], unit=u.deg), None), + (ValueError, cubem, 1*u.s, 1*u.s, 1*u.s), + (ValueError, cubem, u.Quantity([0], unit=u.deg), None, u.Quantity([1.5, 2.], unit=u.deg))]) def test_crop_by_coords_error(test_input): - with pytest.raises(ValueError): - test_input[0].crop_by_coords(*test_input[1:]) + with pytest.raises(test_input[0]): + test_input[1].crop_by_coords(*test_input[2:]) @pytest.mark.parametrize( From 77b32f71388cce25b2f86f7593f034a55bd6d102 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Wed, 25 Apr 2018 13:45:52 +0200 Subject: [PATCH 07/14] New API for crop by coords with units --- ndcube/ndcube.py | 33 ++++++++++++++++++++++++++++----- ndcube/tests/test_ndcube.py | 8 ++++++-- 2 files changed, 34 insertions(+), 7 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index d169b5e2d..b9f3e4df6 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -87,31 +87,37 @@ def world_axis_physical_types(self): pass @abc.abstractmethod - def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): + def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, units=None): """ Crops an NDCube given minimum values and interval widths along axes. Parameters ---------- - lower_corner: iterable of `astropy.units.Quantity` + lower_corner: iterable of `astropy.units.Quantity` or `float` The minimum desired values along each relevant axis after cropping described in physical units consistent with the NDCube's wcs object. The length of the iterable must equal the number of data dimensions and must have the same order as the data. - interval_widths: iterable of `astropy.units.Quantity` + interval_widths: iterable of `astropy.units.Quantity` or `float` The width of the region of interest in each dimension in physical units consistent with the NDCube's wcs object. The length of the iterable must equal the number of data dimensions and must have the same order as the data. This argument will be removed in versions 2.0, please use upper_corner argument. - upper_corner: iterable of `astropy.units.Quantity` + upper_corner: iterable of `astropy.units.Quantity` or `float` The maximum desired values along each relevant axis after cropping described in physical units consistent with the NDCube's wcs object. The length of the iterable must equal the number of data dimensions and must have the same order as the data. + units: iterable of `astropy.units.quantity.Quantity`, optionnal + If the inputs are set without units, the user must set the units + inside this argument as `str`. + The length of the iterable must equal the number of data dimensions + and must have the same order as the data. + Returns ------- result: NDCube @@ -419,7 +425,7 @@ def extra_coords(self): "value": self._extra_coords_wcs_axis[key]["value"]} return result - def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): + def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, units=None): # The docstring is defined in NDDataBase n_dim = self.data.ndim @@ -442,6 +448,23 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None): if (len(lower_corner) != len(upper_corner)) or (len(lower_corner) != n_dim): raise ValueError("lower_corner and upper_corner must have same" "number of elements as number of data dimensions.") + if units: + # Raising a value error if units have not the data dimensions. + if len(units) != n_dim: + raise ValueError('units must have same number of elements as ' + 'number of data dimensions.') + # If inputs are not Quantity objects, they are modified into specified units + for i in range(n_dim): + if type(lower_corner[i]) is not u.quantity.Quantity: + lower_corner[i] = u.Quantity(lower_corner[i], unit=units[i]) + if type(upper_corner[i]) is not u.quantity.Quantity: + upper_corner[i] = u.Quantity(upper_corner[i], unit=units[i]) + else: + for i in range(n_dim): + if type(lower_corner[i]) is not u.quantity.Quantity or \ + type(upper_corner[i]) is not u.quantity.Quantity: + raise ValueError("The inputs are not Quantity objects, you " + "can use the units argument to set it") # Derive all corners coordinates quantity_list = [[lower_corner[i], upper_corner[i]] for i in range(n_dim)] all_corners = [self.world_to_pixel(*a) for a in product(*quantity_list)] diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index 10b004dc2..139cd9f73 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -844,6 +844,8 @@ def test_world_to_pixel(test_input, expected): ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], None, [1.7*u.deg, 1*u.deg, 1.07e-9*u.m]), cubem[:, :, :3]), ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], None, [1*u.s, 2.5*u.arcsec, 0.5*u.arcsec]), + cube_rotated[:, :3, 1:4]), + ((cube_rotated, [0, 1.5, 0], None, [1, 2.5, 0.5], ['s', 'arcsec', 'arcsec']), cube_rotated[:, :3, 1:4])]) def test_crop_by_coords(test_input, expected): helpers.assert_cubes_equal( @@ -852,8 +854,10 @@ def test_crop_by_coords(test_input, expected): @pytest.mark.parametrize("test_input", [ (ValueError, cubem, u.Quantity([0], unit=u.deg), u.Quantity([1.5, 2.], unit=u.deg), None), - (ValueError, cubem, 1*u.s, 1*u.s, 1*u.s), - (ValueError, cubem, u.Quantity([0], unit=u.deg), None, u.Quantity([1.5, 2.], unit=u.deg))]) + (ValueError, cubem, [1*u.s], [1*u.s], [1*u.s]), + (ValueError, cubem, u.Quantity([0], unit=u.deg), None, u.Quantity([1.5, 2.], unit=u.deg)), + (ValueError, cubem, [1], None, [1], ['s', 'deg']), + (ValueError, cubem, [1], None, [1])]) def test_crop_by_coords_error(test_input): with pytest.raises(test_input[0]): test_input[1].crop_by_coords(*test_input[2:]) From 16c08f095a9947290f32fd4a3a77e4e702c5c6b9 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Wed, 25 Apr 2018 14:54:55 +0200 Subject: [PATCH 08/14] Some Quantity types checks refactored --- ndcube/ndcube.py | 14 +++++++------- ndcube/tests/test_ndcube.py | 2 +- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index b9f3e4df6..01bbfe484 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -454,17 +454,17 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, raise ValueError('units must have same number of elements as ' 'number of data dimensions.') # If inputs are not Quantity objects, they are modified into specified units - for i in range(n_dim): - if type(lower_corner[i]) is not u.quantity.Quantity: - lower_corner[i] = u.Quantity(lower_corner[i], unit=units[i]) - if type(upper_corner[i]) is not u.quantity.Quantity: - upper_corner[i] = u.Quantity(upper_corner[i], unit=units[i]) + lower_corner = [u.Quantity(lower_corner[i], unit=units[i]) + for i in range(self.data.ndim)] + upper_corner = [u.Quantity(upper_corner[i], unit=units[i]) + for i in range(self.data.ndim)] else: for i in range(n_dim): if type(lower_corner[i]) is not u.quantity.Quantity or \ type(upper_corner[i]) is not u.quantity.Quantity: - raise ValueError("The inputs are not Quantity objects, you " - "can use the units argument to set it") + raise TypeError("lower_corner and interval_widths/upper_corner must be " + "of type astropy.units.Quantity or the units kwarg " + "must be set.") # Derive all corners coordinates quantity_list = [[lower_corner[i], upper_corner[i]] for i in range(n_dim)] all_corners = [self.world_to_pixel(*a) for a in product(*quantity_list)] diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index 139cd9f73..b61e3a994 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -857,7 +857,7 @@ def test_crop_by_coords(test_input, expected): (ValueError, cubem, [1*u.s], [1*u.s], [1*u.s]), (ValueError, cubem, u.Quantity([0], unit=u.deg), None, u.Quantity([1.5, 2.], unit=u.deg)), (ValueError, cubem, [1], None, [1], ['s', 'deg']), - (ValueError, cubem, [1], None, [1])]) + (TypeError, cubem, [1], None, [1])]) def test_crop_by_coords_error(test_input): with pytest.raises(test_input[0]): test_input[1].crop_by_coords(*test_input[2:]) From 1adde1ddd4bd8cbaa433313a2f721a0d66b4e6a1 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Wed, 25 Apr 2018 15:14:45 +0200 Subject: [PATCH 09/14] Fix a test error --- ndcube/tests/test_ndcube.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index b61e3a994..0977738e2 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -857,7 +857,7 @@ def test_crop_by_coords(test_input, expected): (ValueError, cubem, [1*u.s], [1*u.s], [1*u.s]), (ValueError, cubem, u.Quantity([0], unit=u.deg), None, u.Quantity([1.5, 2.], unit=u.deg)), (ValueError, cubem, [1], None, [1], ['s', 'deg']), - (TypeError, cubem, [1], None, [1])]) + (TypeError, cubem, [1, 2, 3], None, [2, 3, 4])]) def test_crop_by_coords_error(test_input): with pytest.raises(test_input[0]): test_input[1].crop_by_coords(*test_input[2:]) From 76501e13e4db96ba06c1c37365c52413092919d5 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Thu, 26 Apr 2018 11:19:43 +0200 Subject: [PATCH 10/14] Improving units checks in crop by coords --- ndcube/ndcube.py | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index 01bbfe484..75e6ae2f9 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -459,12 +459,10 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, upper_corner = [u.Quantity(upper_corner[i], unit=units[i]) for i in range(self.data.ndim)] else: - for i in range(n_dim): - if type(lower_corner[i]) is not u.quantity.Quantity or \ - type(upper_corner[i]) is not u.quantity.Quantity: - raise TypeError("lower_corner and interval_widths/upper_corner must be " - "of type astropy.units.Quantity or the units kwarg " - "must be set.") + if any([not isinstance(x, u.Quantity) for x in lower_corner + upper_corner]): + raise TypeError("lower_corner and interval_widths/upper_corner must be " + "of type astropy.units.Quantity or the units kwarg " + "must be set.") # Derive all corners coordinates quantity_list = [[lower_corner[i], upper_corner[i]] for i in range(n_dim)] all_corners = [self.world_to_pixel(*a) for a in product(*quantity_list)] From 7c56239444fb069800a0aeccfa00d1bb6f43eebd Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Thu, 26 Apr 2018 16:41:20 +0200 Subject: [PATCH 11/14] Improving corners calculation --- ndcube/ndcube.py | 31 +++++++++++++------------------ 1 file changed, 13 insertions(+), 18 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index 75e6ae2f9..6197b1a6c 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -463,24 +463,19 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, raise TypeError("lower_corner and interval_widths/upper_corner must be " "of type astropy.units.Quantity or the units kwarg " "must be set.") - # Derive all corners coordinates - quantity_list = [[lower_corner[i], upper_corner[i]] for i in range(n_dim)] - all_corners = [self.world_to_pixel(*a) for a in product(*quantity_list)] - # Inputing of all corners coordinates inside a numpy array - # According to the boundary conditions - corners_array = np.zeros((2**n_dim, n_dim)) - for i in range(2**n_dim): - for j in range(n_dim): - corners_array[i, j] = u.Quantity(all_corners[i][j]).value - if corners_array[i, j] > self.data.shape[j]: - corners_array[i, j] = self.data.shape[j] - if corners_array[i, j] < 0: - corners_array[i, j] = 0 - # Taking the maximum and minimum values of coordinates - lower_pixels = corners_array.min(0) - upper_pixels = corners_array.max(0) - # Creating a tuple to crop the data with inputed coordinates - item = tuple([slice(int(lower_pixels[i]), int(upper_pixels[i]) + 1) for i in range(n_dim)]) + # Get all corners of region of interest. + all_world_corners_grid = np.meshgrid(*[u.Quantity([lower_corner[i], upper_corner[i]], + unit=lower_corner[i].unit).value + for i in range(self.data.ndim)]) + all_world_corners = [all_world_corners_grid[i].flatten()*lower_corner[i].unit + for i in range(n_dim)] + # Convert to pixel coordinates + all_pix_corners = self.world_to_pixel(*all_world_corners) + # Derive slicing item with which to slice NDCube. + # Be sure to round down min pixel and round up + 1 the max pixel. + item = tuple([slice(int(axis_pixels.value.min()), + int(np.ceil(axis_pixels.value.max()))+1) + for axis_pixels in all_pix_corners]) return self[item] def crop_by_extra_coord(self, min_coord_value, interval_width, coord_name): From 4f729ddce4d179c0a9911d63db7034283633b525 Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Thu, 26 Apr 2018 16:50:06 +0200 Subject: [PATCH 12/14] Improving the data slicing --- ndcube/ndcube.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index 6197b1a6c..f0943b25d 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -473,7 +473,7 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, all_pix_corners = self.world_to_pixel(*all_world_corners) # Derive slicing item with which to slice NDCube. # Be sure to round down min pixel and round up + 1 the max pixel. - item = tuple([slice(int(axis_pixels.value.min()), + item = tuple([slice(int(np.clip(axis_pixels.value.min(), 0, None)), int(np.ceil(axis_pixels.value.max()))+1) for axis_pixels in all_pix_corners]) return self[item] From b09724169ce4d49f7fe85915cfc9a075397637bd Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Thu, 26 Apr 2018 17:10:57 +0200 Subject: [PATCH 13/14] Fix some error bugs --- ndcube/ndcube.py | 6 +++--- ndcube/tests/test_ndcube.py | 10 +++++----- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index f0943b25d..97f92a701 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -464,9 +464,9 @@ def crop_by_coords(self, lower_corner, interval_widths=None, upper_corner=None, "of type astropy.units.Quantity or the units kwarg " "must be set.") # Get all corners of region of interest. - all_world_corners_grid = np.meshgrid(*[u.Quantity([lower_corner[i], upper_corner[i]], - unit=lower_corner[i].unit).value - for i in range(self.data.ndim)]) + all_world_corners_grid = np.meshgrid( + *[u.Quantity([lower_corner[i], upper_corner[i]], unit=lower_corner[i].unit).value + for i in range(self.data.ndim)]) all_world_corners = [all_world_corners_grid[i].flatten()*lower_corner[i].unit for i in range(n_dim)] # Convert to pixel coordinates diff --git a/ndcube/tests/test_ndcube.py b/ndcube/tests/test_ndcube.py index 0977738e2..8e6e39581 100644 --- a/ndcube/tests/test_ndcube.py +++ b/ndcube/tests/test_ndcube.py @@ -837,16 +837,16 @@ def test_world_to_pixel(test_input, expected): @pytest.mark.parametrize("test_input,expected", [ - ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 5.e-11*u.m], None), + ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], [1*u.deg, 1*u.deg, 4.e-11*u.m], None), cubem[:, :, :3]), ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], [1*u.s, 1*u.arcsec, 0.5*u.arcsec], None), - cube_rotated[:, :3, 1:4]), - ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], None, [1.7*u.deg, 1*u.deg, 1.07e-9*u.m]), + cube_rotated[:, :4, 1:5]), + ((cubem, [0.7*u.deg, 1.3e-5*u.deg, 1.02e-9*u.m], None, [1.7*u.deg, 1*u.deg, 1.06e-9*u.m]), cubem[:, :, :3]), ((cube_rotated, [0*u.s, 1.5*u.arcsec, 0*u.arcsec], None, [1*u.s, 2.5*u.arcsec, 0.5*u.arcsec]), - cube_rotated[:, :3, 1:4]), + cube_rotated[:, :4, 1:5]), ((cube_rotated, [0, 1.5, 0], None, [1, 2.5, 0.5], ['s', 'arcsec', 'arcsec']), - cube_rotated[:, :3, 1:4])]) + cube_rotated[:, :4, 1:5])]) def test_crop_by_coords(test_input, expected): helpers.assert_cubes_equal( test_input[0].crop_by_coords(*test_input[1:]), expected) From d6cf940143ce160420e4e2ccc72f3bd39007f9bc Mon Sep 17 00:00:00 2001 From: Baptiste PELLORCE Date: Fri, 27 Apr 2018 11:10:48 +0200 Subject: [PATCH 14/14] Removing product importation --- ndcube/ndcube.py | 1 - 1 file changed, 1 deletion(-) diff --git a/ndcube/ndcube.py b/ndcube/ndcube.py index 97f92a701..1014898ee 100644 --- a/ndcube/ndcube.py +++ b/ndcube/ndcube.py @@ -2,7 +2,6 @@ import abc import warnings -from itertools import product import numpy as np import astropy.nddata