cta-observatory / ctapipe

Low-level data processing pipeline software for CTAO or similar arrays of Imaging Atmospheric Cherenkov Telescopes
https://ctapipe.readthedocs.org
BSD 3-Clause "New" or "Revised" License
64 stars 268 forks source link

Discrepancy in coordinate transformation in the MuonProcessor #2467

Closed jstvdk closed 10 months ago

jstvdk commented 10 months ago

Describe the bug I think the coordinate transformation in ctapipe.image.muon.MuonProcessor is not properly calculated. To be consistent with the User Guide about Coordinates usage in ctapipe, before performing transformation coordinates to the TelescopeFrame, one should create SkyCoord object in CameraFrame with these coordinates. But in the MuonProcessor we immediately transform the camera.geometry to the TelescopeFrame, and miss the steps of creating the CameraFrame object.

To Reproduce Steps to reproduce the behavior:

  1. Launch MuonProcessor for the muonic event.
  2. Look in the variables fov_lon and fov_lat - they are not consistent with calculations as in the User Guide.

In the current realisation of MuonProcessor they are created in the next way:

geometry = source.subarray.tel[tel_id].camera.geometry.transform_to(TelescopeFrame())
fov_lon = geometry.pix_x
fov_lat = geometry.pix_y

Expected behavior There should be the step of creating the CameraFrame and SkyCoord objects before transformation to the TelescopeFrame, e.g.

camera_frame = CameraFrame(
    focal_length = equivalent_focal_length,
    rotation = source.subarray.tels[tel_id].camera.geometry.cam_rotation,
    telescope_pointing = telescope_pointing,
)
cam_coords = SkyCoord(x = source.subarray.tel[[tel_id].camera.geometry.pix_x,
          y = source.subarray.tel[[tel_id].camera.geometry.pix_x,
          camera_frame,
)

and then we can

telescope_coords = cam_coords.transform_to(TelescopeFrame())
fov_lon = telescope_coords.pix_x
fov_lat = telescope_coords.pix_x

Additional context

This small coordinate discrepancy leads to the difference in optical_efficiency calculations with two methods up to 25%

maxnoe commented 10 months ago

@jstvdk Actually, I just made the opposite fix in #2464

Transforming the CameraGeometry is the correct thing to do, you don't want to essentially re-implement the transformation in the muon code.

I guess that the differences you are seeing is due to effective vs. equivalent focal length.

The code in processor uses the effective focal length by default, which is correct. Due to doing the transformation "by hand" instead of using the transformation of the CameraGeometry, the muon code was always using the equivalent focal length.

So in short, your issue should be fixed by #2464 but in the opposite way you asked here.

jstvdk commented 10 months ago

@maxnoe Thanks for the explanation. I have one more comment, please. I see the changes only in the intensity_fitter.py, but what about processor.py? If I correctly understand, there is a similar transformation in the class MuonProcessor which immediately has an influence on the ring fit itself, because transformed coordinates are used in the ring_fitter method.

maxnoe commented 10 months ago

@jstvdk the fix in #2464 makes it so that the intensity fitter does the same as the processor was already doing.

Basically, the fitter was wrong (due to being written before we introduced CameraGeometry.transform_to and the switch to TelescopeFrame for image parameters), everything else was already correct.

jstvdk commented 10 months ago

@jstvdk the fix in #2464 makes it so that the intensity fitter does the same as the processor was already doing.

Ah, sorry, now I get it. So the problem was because of the different values for focal_length in processor and intensity_fitter, and we should use effective_focal_length in both cases, which is implemented now. Thank you very much for the explanation.