h_max wrt array_level

[thomasgas]

The estimation of h_max in HillasReconstructor: https://github.com/cta-observatory/ctapipe/blob/cc215e2f7233973ff1b604f90c96ac8cea447fcd/ctapipe/reco/HillasReconstructor.py#L333 is giving a result wrt the array observation level. Should we have this number with respect to the sea level? I used the same piece of code used in the hillas_intersector and inside ImPACT to convert the output of the HillasReconstructor (here called reco_hillas) to the "corrected" one:

def correction_xmax(zenith, height, obs_elev):
    height = height*np.cos(zenith)

    # Add on the height of the detector above sea level
    height = height +  obs_elev

    # Lookup this height in the depth tables, the convert Hmax to Xmax
    x_max = thickness_profile(height.to(u.km))
    # Convert to slant depth
    x_max /= np.cos(zenith)
    return x_max

xmax_corrected = correction_xmax(
    np.pi/2 - event.mc.alt.to_value(u.rad), 
    reco_hillas.h_max, 
    event.mcheader.prod_site_alt
)

print(f'reco h_max \t= {reco_hillas.h_max:.2f}')
print(f'mc.energy \t= {event.mc.energy:.2f}')
print(f'from mc \t= {thickness_profile(reco_hillas.h_max.to(u.km)):.2f}')
print(f'event.mc.x_max \t= {event.mc.x_max:.2f}')
print(f'corrected x_max = {xmax_corrected:.2f}')

and this is what i obtain for a huge event, with a good h_max reconstruction (21 telescopes from gamma, Paranal, z20):

reco h_max  = 9035.94 m
mc.energy   = 2.60 TeV
from mc     = 331.01 g / cm2
event.mc.x_max  = 286.86 g / cm2
corrected x_max = 280.40 g / cm2

We correct it or we leave as it is, specifying that it's wrt the array observation level?

[kosack]

I think sea-level makes the most sense (and is probably what anybody would expect it to be). We should change the Field description to say that it is that explicitly.

This also makes me notice that we need to add the site location into the instrument.SubarrayDescription (not just in the MC header, since it needs to be there for real data files as well)

[thomasgas]

This goes together with the additional atmospheric profile for La Palma (https://github.com/cta-observatory/ctapipe-extra/issues/24)