Move get_shower_coordinates to hillas, fixes #934 (#976)

* Move get_shower_coordinates to hillas

ctapipe/image/concentration.py

@@ -1,5 +1,6 @@
 import numpy as np
 from ..io.containers import ConcentrationContainer
+from .hillas import camera_to_shower_coordinates
 
 
 def concentration(geom, image, hillas_parameters):
@@ -21,7 +22,7 @@ def concentration(geom, image, hillas_parameters):
     cog_pixels = np.argsort(delta_x**2 + delta_y**2)
     conc_cog = np.sum(image[cog_pixels[:3]]) / h.intensity
 
-    longi, trans = geom.get_shower_coordinates(h.x, h.y, h.psi)
+    longi, trans = camera_to_shower_coordinates(geom.pix_x, geom.pix_y, h.x, h.y, h.psi)
 
     # get all pixels inside the hillas ellipse
     mask_core = (longi**2 / h.length**2) + (trans**2 / h.width**2) <= 1.0

ctapipe/image/hillas.py

@@ -17,6 +17,43 @@
 ]
 
 
+def camera_to_shower_coordinates(x, y, cog_x, cog_y, psi):
+    '''
+    Return longitudinal and transverse coordinates for x and y
+    for a given set of hillas parameters
+
+    Parameters
+    ----------
+    x: u.Quantity[length]
+        x coordinate in camera coordinates
+    y: u.Quantity[length]
+        y coordinate in camera coordinates
+    cog_x: u.Quantity[length]
+        x coordinate of center of gravity
+    cog_y: u.Quantity[length]
+        y coordinate of center of gravity
+    psi: Angle
+        orientation angle
+
+    Returns
+    -------
+    longitudinal: astropy.units.Quantity
+        longitudinal coordinates (along the shower axis)
+    transverse: astropy.units.Quantity
+        transverse coordinates (perpendicular to the shower axis)
+    '''
+    cos_psi = np.cos(psi)
+    sin_psi = np.sin(psi)
+
+    delta_x = x - cog_x
+    delta_y = y - cog_y
+
+    longi = delta_x * cos_psi + delta_y * sin_psi
+    trans = delta_x * -sin_psi + delta_y * cos_psi
+
+    return longi, trans
+
+
 class HillasParameterizationError(RuntimeError):
     pass
 

ctapipe/image/timing_parameters.py

@@ -4,6 +4,8 @@
 
 import numpy as np
 from ctapipe.io.containers import TimingParametersContainer
+from .hillas import camera_to_shower_coordinates
+
 
 __all__ = [
     'timing_parameters'
@@ -31,26 +33,26 @@ def timing_parameters(geom, image, peakpos, hillas_parameters):
     """
 
     unit = geom.pix_x.unit
-    pix_x = geom.pix_x.value
-    pix_y = geom.pix_y.value
 
     # select only the pixels in the cleaned image that are greater than zero.
     # we need to exclude possible pixels with zero signal after cleaning.
     mask = image > 0
-    pix_x = pix_x[mask]
-    pix_y = pix_y[mask]
+    pix_x = geom.pix_x[mask]
+    pix_y = geom.pix_y[mask]
     image = image[mask]
     peakpos = peakpos[mask]
 
     assert pix_x.shape == image.shape, 'image shape must match geometry'
     assert pix_x.shape == peakpos.shape, 'peakpos shape must match geometry'
 
-    longi, trans = geom.get_shower_coordinates(
+    longi, trans = camera_to_shower_coordinates(
+        pix_x,
+        pix_y,
         hillas_parameters.x,
         hillas_parameters.y,
         hillas_parameters.psi
     )
-    slope, intercept = np.polyfit(longi[mask], peakpos, deg=1, w=np.sqrt(image))
+    slope, intercept = np.polyfit(longi.value, peakpos, deg=1, w=np.sqrt(image))
 
     return TimingParametersContainer(
         slope=slope / unit,

ctapipe/instrument/camera.py

@@ -533,33 +533,6 @@ def get_border_pixel_mask(self, width=1):
         self.border_cache[width] = mask
         return mask
 
-    def get_shower_coordinates(self, x, y, psi):
-        '''
-        Return longitudinal and transverse coordinates of the pixels
-        for a given set of hillas parameters
-
-        Parameters
-        ----------
-        hillas_parameters: ctapipe.io.containers.HilllasContainer
-
-        Returns
-        -------
-        longitudinal: astropy.units.Quantity
-            longitudinal coordinates (along the shower axis)
-        transverse: astropy.units.Quantity
-            transverse coordinates (perpendicular to the shower axis)
-        '''
-        cos_psi = np.cos(psi)
-        sin_psi = np.sin(psi)
-
-        delta_x = self.pix_x - x
-        delta_y = self.pix_y - y
-
-        longi = delta_x * cos_psi + delta_y * sin_psi
-        trans = delta_x * -sin_psi + delta_y * cos_psi
-
-        return longi, trans
-
     def position_to_pix_index(self, x, y):
         '''
         Return the index of a camera pixel which contains a given position (x,y)

docs/tutorials/lst_analysis_bootcamp_2018.ipynb

@@ -465,7 +465,8 @@
    "source": [
     "from ctapipe.image import hillas_parameters, leakage, concentration\n",
     "from ctapipe.image.timing_parameters import timing_parameters\n",
-    "from ctapipe.image.cleaning import number_of_islands"
+    "from ctapipe.image.cleaning import number_of_islands\n",
+    "from ctapipe.image.hillas import camera_to_shower_coordinates"
    ]
   },
   {
@@ -519,7 +520,9 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "long, trans = camera.get_shower_coordinates(hillas.x, hillas.y, hillas.psi)\n",
+    "long, trans = camera_to_shower_coordinates(\n",
+    "    camera.pix_x, camera.pix_y,hillas.x, hillas.y, hillas.psi\n",
+    ")\n",
     "\n",
     "plt.plot(long[clean], dl1.peakpos[0][clean], 'o')\n",
     "plt.plot(long[clean], timing.slope * long[clean] + timing.intercept)"
@@ -1112,7 +1115,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.7.1"
+   "version": "3.6.8"
   }
  },
  "nbformat": 4,