HERA-Team / hera_sim

Simple simulation code for HERA-like redundant interferometric arrays
Other
16 stars 8 forks source link

feat: mutual coupling a la Josaitis+22 #246

Closed r-pascua closed 1 year ago

r-pascua commented 1 year ago

This PR implements the "first order coupling" model described in Josaitis et al. 2022 via the sigchain.MutualCoupling class. Some additional utility functions were also added to help speed things up. This should support fully or partially (linearly) polarized data, although the implementation is a little more memory hungry than it needs to be for unpolarized data (this is a design choice I made arbitrarily—the code is much neater this way).

At a high level, the procedure goes like this:

  1. Calculate a large "coupling matrix"
  2. Reshape the input data array to a shape compatible with the coupling matrix
  3. Calculate the first-order coupling (what I'll call the crosstalk visibilities) via xt_vis = data @ coupling.H + coupling @ vis
  4. Reshape the crosstalk visibilities to comply with the input data array and return the result

This isn't currently ready to be merged; I still need to add documentation and a few unit tests. I'm a bit hesitant to try optimizing further—the only things I can really think of to improve performance require writing some wrapped C-code, or figure out how to access the BLAS/LAPACK functions that numba uses under the hood (assuming I actually have an idea of how numba works...). Here are some specs for timing tests I did for a H4C-sized file with NRAO resources, using 8 cores (I don't think I had any other processes running, but I'm not sure if anything else was running on the same node at the time I ran the tests).

First, the data parameters:

Now for the various steps:

  1. Calculate coupling: ~15 seconds
  2. Reshape input data: ~2.5 seconds
  3. Do matrix products and sum: ~6.5 seconds
  4. Reshape crosstalk visibilities: ~1.3 seconds

Will add some plots at a later date, but preliminary simulated visibilities look good.

steven-murray commented 1 year ago

Thanks @r-pascua, this is really impressive. With the timing, I suppose the killer will be the matrix product (I assume the initial 15 sec for computing coupling matrix is a one-off for the array, and that each integration uses the same coupling matrix?). Multiplying by ~8000 is only about 14 hours, which is very doable. At that point it is almost certainly not the bottleneck any more.

codecov[bot] commented 1 year ago

Codecov Report

Base: 96.54% // Head: 96.52% // Decreases project coverage by -0.03% :warning:

Coverage data is based on head (ac9fed1) compared to base (3a2c1f8). Patch coverage: 96.71% of modified lines in pull request are covered.

Additional details and impacted files ```diff @@ Coverage Diff @@ ## main #246 +/- ## ========================================== - Coverage 96.54% 96.52% -0.03% ========================================== Files 23 24 +1 Lines 2808 3019 +211 ========================================== + Hits 2711 2914 +203 - Misses 97 105 +8 ``` | [Impacted Files](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team) | Coverage Δ | | |---|---|---| | [hera\_sim/utils.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vdXRpbHMucHk=) | `97.07% <93.93%> (-2.07%)` | :arrow_down: | | [hera\_sim/simulate.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vc2ltdWxhdGUucHk=) | `94.72% <95.65%> (-0.25%)` | :arrow_down: | | [hera\_sim/sigchain.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vc2lnY2hhaW4ucHk=) | `98.37% <98.37%> (-0.03%)` | :arrow_down: | | [hera\_sim/components.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vY29tcG9uZW50cy5weQ==) | `91.91% <100.00%> (+0.08%)` | :arrow_up: | | [hera\_sim/eor.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vZW9yLnB5) | `100.00% <100.00%> (ø)` | | | [hera\_sim/foregrounds.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vZm9yZWdyb3VuZHMucHk=) | `100.00% <100.00%> (ø)` | | | [hera\_sim/interpolators.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vaW50ZXJwb2xhdG9ycy5weQ==) | `97.05% <100.00%> (ø)` | | | [hera\_sim/noise.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vbm9pc2UucHk=) | `100.00% <100.00%> (ø)` | | | [hera\_sim/rfi.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vcmZpLnB5) | `100.00% <100.00%> (ø)` | | | [hera\_sim/visibilities/simulators.py](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team#diff-aGVyYV9zaW0vdmlzaWJpbGl0aWVzL3NpbXVsYXRvcnMucHk=) | `85.46% <100.00%> (ø)` | | | ... and [1 more](https://codecov.io/gh/HERA-Team/hera_sim/pull/246?src=pr&el=tree-more&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team) | | Help us with your feedback. Take ten seconds to tell us [how you rate us](https://about.codecov.io/nps?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team). Have a feature suggestion? [Share it here.](https://app.codecov.io/gh/feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=HERA-Team)

:umbrella: View full report at Codecov.
:loudspeaker: Do you have feedback about the report comment? Let us know in this issue.

review-notebook-app[bot] commented 1 year ago

Check out this pull request on  ReviewNB

See visual diffs & provide feedback on Jupyter Notebooks.


Powered by ReviewNB

steven-murray commented 1 year ago

Thanks @r-pascua , this is excellent work. I would love to see a notebook demo with some plots that "show" what the tests test. If we could go through such a notebook on a validation call, I think that would be the best final review. Otherwise, the code changes look good to me.

Funnily enough, I thought I saw such a notebook in the files changed in this PR, but now I can't see it...

r-pascua commented 1 year ago

Thanks for the review! I do have a notebook in the works, but it was a little messy and I wasn't sure if it was necessary for this PR. I'll have it ready by next week's validation call.