Align Einstein Telescope geometry with LAL by milanwils · Pull Request #894 · bilby-dev/bilby

milanwils · 2025-01-23T09:43:47Z

Problem
As mentioned in #713, the current geometry of the ET is not in line with its definition in LAL.

The labels of ET2 and ET3 are swapped.
Currently the new corner points are calculated by updating the latitude and longitude approximately and not updating the elevation. Additionally, the azimuth and tilt do not take into account the curvature of the earth. This leads to small differences in the vertex position and unit vector. As a result, the null stream is not accurate and signal leaks into it.

Proposed Solution

This is addressed in this PR by adding an optional clockwise option to the detector config. I set the value in the ET.interferometer file to clockwise=False so that the default is consistent with LAL. If this is undesired for backwards compatibility, it could be set to clockwise=True .
This issue is addressed by following the same specification as LAL: T1400308. The code is rather verbose to avoid confusion but can be changed if denser code is preferred. This algorithm required the implementation of a transform from WGS84 coordinated to ECEF coordinates. I added a small test to check the consistency with the forward transformation. In the implementation of this transform I used uninformative variable names to be consistent with the reference (they do not have a physical meaning anyway).

Additional information
There is also a small error in the LAL config, see LAL PR 2399 . Only if both get accepted do the configurations become identical.

ColmTalbot

This is a nice change, I think we'll need to think a little about the best way to do this. I've suggested waiting a while to get this in as I don't think it is urgent, the labelling of ET2/ET3 doesn't really matter, and we can make the changes more satisfyingly if we break things.

ColmTalbot · 2025-01-23T16:07:57Z

-                * 180
-                / np.pi
-            )
+            unit_vector_x = self[ii].geometry.unit_vector_along_arm("x")


I know this is not the change for the current code, but given that this is clearly better thought out than the current method, can the logic for updating the geometry be moved into a standalone function (with a nice docstring)?

Perhaps we could add setters to the InterferometerGeometry class for x, y and vertex? For the lat, long, elevation, tilt and azimuth these setters are separate but it does not seem sensible to me that someone can change only one Cartesian coordinate.

All the update logic could then be moved into the geometry class and only the calculation of the new corner points (which is triangle specific) remains in the 'TriangularInterferometer` class.

My main concern with this is approach is that one would need to be careful to update the vertex position before changing the unit vector. Otherwise the calculation of tilt and azimuth is wrong. This could be avoided by making the unit vectors the "base information" and deriving the tilt and azimuth from it.

This also raises the question of whether the constructor should be adapted to allow users to define the detector directly from Cartesian parameters. I guess you could argue that this is actually the only way you should be able to set these parameters as interferometers are constant.

ColmTalbot · 2025-01-23T16:10:03Z

        yarm_azimuth,
        xarm_tilt=0.0,
        yarm_tilt=0.0,
+        clockwise=True,


I would personally vote for just making a breaking change here (and waiting for a major release) to avoid adding a new parameter that doesn't really mean anything (the set of three Interferometers are just a permutation of each other with the different directions).

Fine by me, if people really want they can still override the detector names

ColmTalbot · 2025-01-23T16:15:16Z

    return np.array([x_comp, y_comp, z_comp])


+def get_vertex_position_ellipsoid(x_comp, y_comp, z_comp):


I feel like this function should take and recieve the fame type, either (float, float, float) -> (float, float, float) or array -> array. I don't really mind which, but it would be good to be consistent.

I now see that the same thing happens in get_vertex_position_geocentric. sigh, I'm not sure if it is worth changing the other function.

We could change the other function to something like

def get_vertex_position_geocentric(position, _longitude=None, _elevation=None): if longitude is not None and not isarray(position): WARN("... syntax changed in Bilby 3.0.0 ...") position = np.array([position, _longitude, _elevation]) ...

This would need a wider set of eyes/opinions.

I went through the exact same thought process when writing the function. This seems like a good compromise.

milanwils · 2025-03-06T10:14:34Z

Update: LAL MR2399 has been approved and the methodology in T1400308 now contains the algorithm added in this PR.

milanwils · 2025-10-30T14:17:20Z

For future reference: the transformation is also implemented in astropy (this is what pycbc uses).

mj-will · 2026-08-26T08:35:48Z

Now we're preparing for 3.0, I think we can more easily consider the breaking changes mentioned above. @milanwils @ColmTalbot can you revisit this and assess where we stand?

milanwils · 2026-09-14T10:56:44Z

@mj-will Sorry for the delay in replying. I had to refresh my memory a bit because I made this request quite some time ago. As far as I can see, the fix in #1016 solved most of the issues but to be entirely correct, the changes in this PR should still be applied because:

The distance parameter is used as the distance along the surface of the earth, not the distance of the detector arm
The elevation is ignored in the calculations
The azimuth is only updated approximately
The tilt is not updated (local plane at each corner point is slightly different)

I believe there are just a few questions to be answered and then we can proceed rather quickly:

Do we want to make the labelling of the detectors consistent with LAL? I believe this is the only breaking change.
Do we remove the clockwise argument I introduced? I suggest yes, given that @ColmTalbot was in favour.
What do we do with the function signatures of e.g. get_vertex_position_ellipsoid? For me it is fine to use @ColmTalbot 's suggestion.
Do we want to use astropy for the conversion of coordinates or the implementation I made?
Where to put the logic for calculating the next detector parameters. Currently this is in the initialiser but this can be put into a member function of the class. This would return 7 parameters, which might be a bit cumbersome. I would personally just prefer a decent docstring with the initialiser function.

mj-will · 2026-09-14T12:58:39Z

Thanks @milanwils, here are my thoughts.

1. Do we want to make the labelling of the detectors consistent with LAL? I believe this is the only breaking change.

Yes, I think so.

2. Do we remove the clockwise argument I introduced? I suggest yes, given that @ColmTalbot  was in favour.

Agree, this would be a breaking change but that's fine for 3.0

3. What do we do with the function signatures of e.g. `get_vertex_position_ellipsoid`? For me it is fine to use @ColmTalbot 's suggestion.

Let's start with that and see

4. Do we want to use astropy for the conversion of coordinates or the implementation I made?

My concern with astropy would be speed. Do you know how they compare? If your implementation is faster then a test to show they're equivalent may be enough.

5. Where to put the logic for calculating the next detector parameters. Currently this is in the initialiser but this can be put into a member function of the class. This would return 7 parameters, which might be a bit cumbersome. I would personally just prefer a decent docstring with the initialiser function.

I've not looked at the PR in detail, but if I'm following correctly, you're referring to the large block of code that's been added to the the initialiser. My preference would be to try and factor this out into a method (or methods) so it's easier to test in isolation but this might need some more thought.

milanwils · 2026-09-16T11:15:46Z

@mj-will

I have pushed changes for 1-3.

Is speed really important here? This function is only supposed to be called when creating the detector object or when debugging I suppose. Even if it would be a fraction slower, I don't see how that would cause any problem.
Ok, I will try to move this logic somewhere else.

ColmTalbot · 2026-09-16T11:17:51Z

I agree with all Michael's comments.
Especially moving the new code in the initializer into a new method, or even a standalone function in gw/utils or gw/geometry that has been added since this PR opened.

The other reason not to use astropy is that they don't yet support array backend flexibility. We should add that to this change. It should be straightforward as the scipy functionality we use is supported, so it should just be sprinkling to xp.

We have docs here and I'm happy to give specific pointers if that's better, it will just take a little while.

milanwils · 2026-09-16T16:15:36Z

I took a look at the array backend flexibility and it should indeed not be hard to implement. To test that everything is working I have to either rebase onto or merge main into this feature branch because the branches diverted before backend flexibility was added to Bilby. What is the standard on this project?

ColmTalbot · 2026-09-16T16:23:18Z

I usually try to rebase, but if it is painful, I just merge. We usually squash commits anyway before we merge anyway, so it shouldn't really matter.

milanwils · 2026-09-17T13:36:16Z

I went with a rebase because it was quite painful.

I stuck with my own implementation for the coordinate transform, moved the update logic to a member function of the TriangularInterferometer class (the update is specific to such a design) and I have added some array back-end flexibility code that works to the best of my knowledge. For the operations on scalars, I left in the numpy code.

As far as I can see, all requested changes are now included. Let me know if there are more changes I have to make. I will still verify that LAL and Bilby are now in perfect alignment.

Are there additional tests we want to add, e.g. is the sum of the antenna patterns zero?

milanwils · 2026-09-17T16:04:09Z

test_et_lal_consistency.py
The new implementation is identical to the LAL implementation as expected.

ColmTalbot

Thanks! Just a few comments about the array backend stuff.

ColmTalbot · 2026-09-21T12:50:32Z

+    semi_major_axis = 6378137  # for ellipsoid model of Earth, in m
+    semi_minor_axis = 6356752.314  # in m
+
+    r = np.sqrt(position[0] ** 2 + position[1] ** 2)


All of the uses of np in this function (and the others) should be xp.

Included in the last commit

ColmTalbot · 2026-09-21T12:53:03Z

+        next_vertex_ellipsoid = get_vertex_position_ellipsoid(next_vertex_geocentric)
+        next_latitude_rad, next_longitude_rad, next_elevation = next_vertex_ellipsoid
+
+        rotation_vector = np.cross(unit_vector_x, unit_vector_y)


As with other cases, these should all be xp.

Included in the last commit

ColmTalbot · 2026-09-21T12:56:06Z

        if not self._vertex_updated:
-            self._vertex = gwutils.get_vertex_position_geocentric(self._latitude, self._longitude,
-                                                                  self.elevation)
+            self._vertex = gwutils.get_vertex_position_geocentric(np.array([self._latitude, self._longitude,


I'm not sure the best way to find the appropriate array backend here, but I would guess we want this import at the top of the file

from ...compat.utils import array_module

then here we would have

xp = array_module((self._latitude, self._longitude, self.elevation))

I have decided to take the approach mentioned in the commit message:

To infer the correct array backend in InterferometerGeometry.vertex,
the property array_backend was added in line with the implementations
in Interferometer and InterferometerList. This also makes the
Interferometer.array_backend method independent of implementation
details in InterferometerGeometry.

I believe this is the cleanest solution, but if you prefer something else then that is fine by me too.

Unrelated question: is there is reason the methods like set_array_backend are used rather than a decorated

@array_backend.setter def array_backend(self, xp): <implementation details>

like there is in many other parts of the code?

No, I guess just a moment of madness on my side, we should probably change this to be a setter as you suggest.

To infer the correct array backend in `InterferometerGeometry.vertex`, the property `array_backend` was added in line with the implementations in `Interferometer` and `InterferometerList`. This also makes the `Interferometer.array_backend` method independent of implementation details in `InterferometerGeometry`.

ColmTalbot

Thanks, just one last place where there is an np instead of xp.

ColmTalbot · 2026-09-30T14:03:36Z

        if not self._vertex_updated:
-            self._vertex = gwutils.get_vertex_position_geocentric(self._latitude, self._longitude,
-                                                                  self.elevation)
+            self._vertex = gwutils.get_vertex_position_geocentric(np.array([self._latitude, self._longitude,


No, I guess just a moment of madness on my side, we should probably change this to be a setter as you suggest.

Co-authored-by: Colm Talbot <talbotcolm@gmail.com>

…r property

milanwils · 2026-09-30T15:42:20Z

Sorry for missing that np call. This time I checked that at least I did not introduce any new calls to np in this PR by running git diff main | grep -C20 "\+.*np".

I also changed the setter because it was not to much work. If you prefer that in a separate PR then I'll revert the commit and open a new PR.

ColmTalbot

Thanks @milanwils! As long as the tests pass, I'm happy with this now.

ColmTalbot · 2026-10-02T15:34:36Z

        if not self._vertex_updated:
-            self._vertex = gwutils.get_vertex_position_geocentric(self._latitude, self._longitude,
-                                                                  self.elevation)
+            self._vertex = gwutils.get_vertex_position_geocentric(xp.array([self._latitude, self._longitude,


I should have caught these, they should all be xp.asarray to appease torch.

milanwils · 2026-10-05T08:27:12Z

@ColmTalbot Sorry, I missed the line in the contributing guideline about testing with multiple backends, I thought this was done automatically. I pushed the changes .

One test was failing in `gw/utils.py` because `torch.deg2rad` only accepts Tensor and not a float. This is fixed by turning the argument into an array using `xp.asarray`. Additionally, the test used assertAlmostEqual which requires the builtin function `__round__` which is not defined for `Tensor` objects. Therefore, the tensors are cast into float before comparing.

milanwils · 2026-10-06T11:58:05Z

@ColmTalbot I made some changes in the test so that it also works for pytorch which forces Tensor inputs (no float) and does not have the builtin round.

ColmTalbot

LGTM, thanks for pushing this @milanwils

Since bilby-dev/bilby#1016 each next channel sits at the end of the current channel's x arm, but the arm azimuths are still rotated by 240 degrees, so the second and third channels have one arm 60 degrees off the triangle (and swapped labels with respect to LAL). Follow LAL (LIGO-T1400308) and bilby-dev/bilby#894 instead: step along the straight x arm on the WGS84 ellipsoid and rotate both arms by 120 degrees about the detector normal, with azimuths, tilts and elevation recomputed in the next vertex's local frame. The colocated option rotates the arms in the same way about the first vertex. The geocentric-to-geodetic inversion uses |z| and real cube roots, so it also holds in the southern hemisphere, where #894's version returns NaN. For the ET-EMR v70 injection this swaps the ET-EMR2/3 SNRs and changes the weaker one from 162.30 to 160.84. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>

milan-wils added 4 commits January 16, 2025 14:25

naming ET detectors anti-clockwise

a5f498e

Implemented transform ECEF to WGS84

ca46c1e

added test for ECEF to WGS84 conversion

3fa49f0

Changed triangular config in line with LAL definition of ET

058d588

ColmTalbot reviewed Jan 23, 2025

View reviewed changes

ColmTalbot mentioned this pull request Nov 5, 2025

BUG: fix triangular interferometer coordinates #1016

Merged

ColmTalbot added bug Something isn't working detector labels Nov 6, 2025

mj-will added waiting for author input waiting on review labels Aug 26, 2026

mj-will linked an issue Aug 26, 2026 that may be closed by this pull request

Einstein Telescope's geometry in Bilby is different from definitions in LAL #713

Open

detector: removed the clockwise option in TriangularInterferometer

c12e576

Changed function signature of get_vertex_* functions

97cc3cd

milanwils force-pushed the et-config branch from 3acb12f to 97cc3cd Compare September 16, 2026 16:00

milan-wils added 2 commits September 16, 2026 18:39

Merge branch 'main' into et-config

fac62ca

moved the update logic to a separate function

7d65481

milan-wils added 2 commits September 17, 2026 16:36

style: pre-commit checks

e7bc7a6

fix: syntax error for line continuation

3ab8e06

ColmTalbot requested changes Sep 21, 2026

View reviewed changes

ColmTalbot reviewed Sep 30, 2026

View reviewed changes

milanwils and others added 3 commits September 30, 2026 17:04

Apply suggestion from @ColmTalbot about remaining np function

60b2c46

Co-authored-by: Colm Talbot <talbotcolm@gmail.com>

changed two remaining np calls to xp

6847347

refactor: replace set_array_backend methods with @array_backend.sette…

2653318

…r property

ColmTalbot approved these changes Oct 2, 2026

View reviewed changes

ColmTalbot reviewed Oct 2, 2026

View reviewed changes

ColmTalbot requested changes Oct 2, 2026

View reviewed changes

fix: changed xp.array to xp.asarray for torch compatibility

b72bc0f

ColmTalbot approved these changes Oct 7, 2026

View reviewed changes

ColmTalbot requested a review from a team October 7, 2026 17:40

		return np.array([x_comp, y_comp, z_comp])


		def get_vertex_position_ellipsoid(x_comp, y_comp, z_comp):

Conversation

milanwils commented Jan 23, 2025

Uh oh!

ColmTalbot left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

milanwils commented Mar 6, 2025

Uh oh!

milanwils commented Oct 30, 2025

Uh oh!

mj-will commented Aug 26, 2026

Uh oh!

milanwils commented Sep 14, 2026

Uh oh!

mj-will commented Sep 14, 2026

Uh oh!

milanwils commented Sep 16, 2026

Uh oh!

ColmTalbot commented Sep 16, 2026

Uh oh!

milanwils commented Sep 16, 2026

Uh oh!

ColmTalbot commented Sep 16, 2026

Uh oh!

milanwils commented Sep 17, 2026

Uh oh!

milanwils commented Sep 17, 2026

Uh oh!

ColmTalbot left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

ColmTalbot left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

milanwils commented Sep 30, 2026

Uh oh!

ColmTalbot left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

milanwils commented Oct 5, 2026

Uh oh!

milanwils commented Oct 6, 2026

Uh oh!

ColmTalbot left a comment

Choose a reason for hiding this comment

Uh oh!

This branch has not been deployed

Reviewers