Skip to content

Commit

Permalink
Few bugs corrected in HillasIntersection after #896
Browse files Browse the repository at this point in the history
  • Loading branch information
Thomas Gasparetto committed Feb 1, 2019
1 parent 37f70ce commit 8a9b821
Showing 1 changed file with 14 additions and 12 deletions.
26 changes: 14 additions & 12 deletions ctapipe/reco/hillas_intersection.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,18 +76,22 @@ def predict(self, hillas_parameters, tel_x, tel_y, array_direction):
src_x, src_y, err_x, err_y = self.reconstruct_nominal(hillas_parameters)
core_x, core_y, core_err_x, core_err_y = self.reconstruct_tilted(
hillas_parameters, tel_x, tel_y)
err_x *= u.rad
err_y *= u.rad
err_x *= u.deg
err_y *= u.deg

nom = SkyCoord(
x=src_x * u.rad,
y=src_y * u.rad,
frame=NominalFrame(array_direction=array_direction)
sky_pos = SkyCoord(
az=src_x * u.deg,
alt=src_y * u.deg,
frame=HorizonFrame()
)
horiz = nom.transform_to(HorizonFrame())

nom_frame = NominalFrame(origin=array_direction)

nom = sky_pos.transform_to(nom_frame)

result = ReconstructedShowerContainer()
result.alt, result.az = horiz.alt, horiz.az
result.alt = nom.altaz.alt # src_y * u.deg

This comment has been minimized.

Sign in to view

Copy link
@maxnoe

maxnoe Feb 2, 2019 •
edited
Loading

Member

What is nom.altaz?

This comment has been minimized.

Sign in to view

Copy link
@thomasgas

thomasgas Feb 3, 2019

If i'm not wrong, the result of self.reconstruct_nominal(hillas_parameters) is in degrees in the HorizonFrame and here I just put this result in the HorizonFrame and transform to the NominalFrame. I could avoid this step just doing:

result.alt = src_y * u.deg
result.az = src_x * u.deg

This comment has been minimized.

Sign in to view

Copy link
@maxnoe

maxnoe Feb 3, 2019

Member

No, I meant what is the altatz member of a NominalFrame coordinate? What does it do, mean?

This comment has been minimized.

Sign in to view

Copy link
@thomasgas

thomasgas Feb 3, 2019

nom gives the coordinates back in the NominalFrame (and you can recover them with nom.delta_az and nom.delta_alt), but with nom.altaz I get the coordinates in the HorizonFrame. Is this what you are asking?

This comment has been minimized.

Sign in to view

Copy link
@maxnoe

maxnoe Feb 4, 2019

Member

I'm asking why and if this is really the case, because I didn't intentionally implement it like this. So where does this feature come from and is it really the nom coordinate transformed to altaz?

This comment has been minimized.

Sign in to view

Copy link
@kosack

kosack Feb 4, 2019

Contributor

that comes from SkyCoord, which has quick helper attributes for common transforms. The .altaz attribute is the same as calling coord.transform_to(c.AltAz). I'm not sure where it assumes the time and location come from though - worth verifying that it does what we expect.

This comment has been minimized.

Sign in to view

Copy link
@maxnoe

maxnoe Feb 4, 2019

Member

It's only chance is to take it from the SkyCoord instance from which it is being called.

This comment has been minimized.

Sign in to view

Copy link
@thomasgas

thomasgas Feb 4, 2019

so in this code here, array_direction is a HorizonFrame with obstime and location specified. Then sky_pos has no info on obstime and location because it's given an HorizonFrame(), with no parameters.
Then nom_frame is this:

<NominalFrame Frame (origin=<AltAz Coordinate (obstime=2018-11-01T02:00:00.000, 
location=(5327448.995782901, -1718665.7386956865, 3051566.9029540345) m, pressure=0.0 hPa, 
temperature=0.0 deg_C, relative_humidity=0, obswl=1.0 micron): (az, alt) in deg
    (180.00000501, 69.99999967)>, obstime=None, location=None)

and nom is this

<SkyCoord (NominalFrame: origin=<AltAz Coordinate (obstime=2018-11-01T02:00:00.000, 
location=(5327448.995782901, -1718665.7386956865, 3051566.9029540345) m, pressure=0.0 hPa, 
temperature=0.0 deg_C, relative_humidity=0, obswl=1.0 micron): (az, alt) in deg
    (180.00000501, 69.99999967)>, obstime=None, location=None): (delta_az, delta_alt) in deg
    (-0.03219597, -0.03347852)>

So the obstime and location are not correctly passed from one system to the other one.

--> One way to correct this behaviour might be to create sky_pos in a HorizonFrame defined in this way:

HorizonFrame(obstime=array_direction.obstime, location=array_direction.location)

So that nom becomes:

<SkyCoord (NominalFrame: origin=<AltAz Coordinate (obstime=2018-11-01T02:00:00.000, 
location=(5327448.995782901, -1718665.7386956865, 3051566.9029540345) m, pressure=0.0 hPa, 
temperature=0.0 deg_C, relative_humidity=0, obswl=1.0 micron): (az, alt) in deg
    (180.00000501, 69.99999967)>, obstime=2018-11-01T02:00:00.000, location=(5327448.995782901, -1718665.7386956865, 3051566.9029540345) m): (delta_az, delta_alt) in deg
    (-0.03219597, -0.03347852)>

This comment has been minimized.

Sign in to view

Copy link
@maxnoe

maxnoe Feb 4, 2019

Member

Yes, create the horizon frame with obstime and location.

Also: I think we should just drop HorizonFrame and replace it by just using AltAz directly from astropy.coordinates.
result.az = nom.altaz.az # src_x * u.deg

tilt = SkyCoord(
x=core_x * u.m,
Expand All @@ -99,7 +103,7 @@ def predict(self, hillas_parameters, tel_x, tel_y, array_direction):
result.core_y = grd.y

x_max = self.reconstruct_xmax(
nom.x, nom.y,
nom.altaz.az, nom.altaz.alt,
tilt.x, tilt.y,
hillas_parameters,
tel_x, tel_y,
Expand Down Expand Up @@ -135,7 +139,7 @@ def reconstruct_nominal(self, hillas_parameters, weighting="Konrad"):
Returns
-------
Reconstructed event position in the nominal system
Reconstructed event position in the horizon system
"""
if len(hillas_parameters) < 2:
Expand Down Expand Up @@ -185,8 +189,6 @@ def reconstruct_nominal(self, hillas_parameters, weighting="Konrad"):
var_x = np.average((sx - x_pos) ** 2, weights=weight)
var_y = np.average((sy - y_pos) ** 2, weights=weight)

# Copy into nominal coordinate

return x_pos, y_pos, np.sqrt(var_x), np.sqrt(var_y)

def reconstruct_tilted(self, hillas_parameters, tel_x, tel_y,
Expand Down

0 comments on commit 8a9b821

Please sign in to comment.