CameraGeometry position_to_pix_index exceptions handling

[mexanick]

position_to_pix_index method handles identically cases where the position is outside of the camera and cases where position is in the gap between pixels. In both cases warning Coordinate (x,y) lies outside camera is raised.

Can this be changed to differentiate the outlined cases?

[maxnoe]

Maybe just a better error message? No pixel at position ...? Or do you really need to be able to distinguish these two cases?

It's not that easy for the irregular camera geometries we have

[kosack]

I guess if it's really required, we could construct a convex hull using the vertices of the boundary pixels of each camera and test if we are inside or outside, but I agree with @maxnoe that maybe it's not really necessary? Or is there a strong use case for knowing if you are in a gap vs outside the camera?

[maxnoe]

This could use shapely as demonstrated here:

https://github.com/cta-observatory/ctapipe/issues/855#issuecomment-496578056

[maxnoe]

import matplotlib.pyplot as plt
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection
from matplotlib.path import Path
import numpy as np
from shapely.geometry import MultiPolygon, Polygon as ShapelyPolygon
from shapely.ops import unary_union

from ctapipe.instrument import CameraGeometry
from ctapipe.visualization import CameraDisplay

cameras = ["LSTCam", "CHEC"]

dummy_fig, dummy_ax = plt.subplots(1, 1)

for area_scale in [1.0, 1.2, 1.8]:
    fig, axs = plt.subplots(1, len(cameras), figsize=(len(cameras) * 5, 5), constrained_layout=True)

    for ax, camera_name in zip(axs.ravel(), cameras):
        camera = CameraGeometry.from_name(camera_name)

        camera.pix_area *= area_scale
        geom = CameraDisplay(camera, ax=dummy_ax)
        pixels = geom.pixels

        patches = []

        polygons = [ShapelyPolygon(p.vertices) for p in pixels.get_paths()]
        boundary = unary_union(polygons)

        if isinstance(boundary, MultiPolygon):
            patches = [
                Polygon(np.array(polygon.exterior.coords))
                for polygon in boundary.geoms
            ]
        else:
            patches = [Polygon(boundary.exterior.coords)]

        col = PatchCollection(patches)
        col.set_cmap('inferno')
        col.set_array(np.arange(len(patches)))
        col.set_clim(0, len(patches))
        col.set_edgecolor('C0')
        col.set_linewidth(3)

        ax.add_collection(col)
        ax.margins(0.05)
        ax.set_aspect(1)
        ax.set_title(camera_name)

    fig.suptitle(f'pixel area scaled up by: {area_scale}')
    fig.savefig(f'whole_camera_{area_scale}.png', dpi=200)
    plt.close(dummy_fig)
plt.show()

[kosack]

You could use scipy.spatial instead to avoid a new dependency on shapely. Requires a bit more code, but doesn't seem too hard. I found this example on StackOverflow:

https://stackoverflow.com/questions/16750618/whats-an-efficient-way-to-find-if-a-point-lies-in-the-convex-hull-of-a-point-cl

[maxnoe]

But we don't want the convex hull, right? Just the hull. Which requires knowledge about the connections, not just the points.

[maxnoe]

I.e. the convex hull for chec would include a diagonal line from the outer edges of the modules like this:

[mexanick]

For my particular use case a convex hull based differentiation would be enough. And differentiating between "in (convex) hull" vs "out of (convex) hull" is desirable, but not too much important for me. As of now I simply ignore the raised warnings, however it is a sort of unsafe practice, as I know for a fact that in my case the normal behavior is still to have a point inside the hull.