Repository: astropy/astropy

Issue:
FITS HDUList deepcopy / BinTableHDU.copy() appears to leak memory
### Description

Doing `deepcopy()` fits `HDUList` appears to leak memory. See the test script below.

The motivational scenario is https://github.com/lightkurve/lightkurve/issues/1388 . `lightkurve` v2.4.2 introduced a workaround involving `deepcopy(hdulist)`, which appears to leak memory.


### Expected behavior

Memory usage should remain stable.

### How to Reproduce

Run the test script below, which
1. downloads a FITS file 
  a. it's a TESS lightcurve FITS file. A variety of FITS file from TESS / Kepler mission have been tried with the same behavior.
3. repeatedly opens the FITS file as HDUList, do a `deepcopy()`, and measure memory usage
  a. the `HDUList` opened is properly closed.
  b. the real life scenario is batch processing a large number of FITS files. The behavior is the same as repeatedly processing the same one in the test script.
4. shows the memory usage keeps increasing. 

```python
from astropy.utils.data import download_file
from astropy.io import fits

from copy import deepcopy
import gc
import tracemalloc

# One TESS lightcurve FITS file from
# https://archive.stsci.edu/missions/tess/download_scripts/sector/tesscurl_sector_47_lc.sh
fits_url = (
    "https://mast.stsci.edu/api/v0.1/Download/file/?uri="
    "mast:TESS/product/tess2021364111932-s0047-0000000014431677-0218-s_lc.fits"
)
iterations_warmup = 40
iterations_run = 60


def read_lite(filename):
    # simulate the logic in
    # https://github.com/lightkurve/lightkurve/pull/1299
    with fits.open(filename) as hdulist:
        # hdulist = hdulist.copy()  # no mem leak
        hdulist = deepcopy(hdulist)  # show mem leak
    return hdulist


print("Downloading fits from: ", fits_url)
fits_local_path = download_file(fits_url, cache=True)
print("Loading fits from: ", fits_local_path)

tracemalloc.start()
print("Warming up...")
for i in range(iterations_warmup):
    hdul = read_lite(fits_local_path)
    hdul.close()
    # print memory used so far in MB
    cur, peak = tracemalloc.get_traced_memory()
    print(i, cur / 1e6, peak / 1e6)

gc.collect()
print("\nAfter GC post warm up:")
cur, peak = tracemalloc.get_traced_memory()
print("  ", cur / 1e6, peak / 1e6)


print("\nRun...")
for i in range(iterations_warmup, iterations_warmup + iterations_run):
    hdul = read_lite(fits_local_path)
    hdul.close()
    # print memory used so far in MB
    cur, peak = tracemalloc.get_traced_memory()
    print(i, cur / 1e6, peak / 1e6)

gc.collect()
print("\nAfter GC post run:")
cur, peak = tracemalloc.get_traced_memory()
print("  ", cur / 1e6, peak / 1e6)
```

Output:
```
Warming up...
0 15.452805 25.355049
1 17.531876 25.355049
2 19.596948 25.355049
...
37 91.717726 93.983287
38 93.780519 96.04608
39 95.843666 98.108921
After GC post warm up:
   95.779748 98.108921

Run...
40 97.903508 100.163626
41 99.966449 102.232118
42 102.03016 104.295829
...
97 215.47148 217.737095
98 217.534155 219.799716
99 219.59678 221.862395
After GC post run:
   219.530556 221.862395
```


### Versions

- Tried on Astropy 5.3.4 and latest `main` at 8b9ca69e7c .
- Unlikely to be OS dependent: the motivational scenario happened on Linux; reproduced with the narrow case here on Windows 11.

```
Windows-10-10.0.22621-SP0
Python 3.10.13 | packaged by conda-forge | (main, Oct 26 2023, 18:01:37) [MSC v.1935 64 bit (AMD64)]
astropy 6.1.dev117+g8b9ca69e7c
Numpy 1.26.2
pyerfa 2.0.1.1
```

Discussion:
github-actions[bot]: Thank you for your contribution to Astropy! 🌌 This checklist is meant to remind the package maintainers who will review this pull request of some common things to look for.

  - [ ] Do the proposed changes actually accomplish desired goals?
  - [ ] Do the proposed changes follow the [Astropy coding guidelines](https://docs.astropy.org/en/latest/development/codeguide.html)?
  - [ ] Are tests added/updated as required? If so, do they follow the [Astropy testing guidelines](https://docs.astropy.org/en/latest/development/testguide.html)?
  - [ ] Are docs added/updated as required? If so, do they follow the [Astropy documentation guidelines](https://docs.astropy.org/en/latest/development/docguide.html#astropy-documentation-rules-and-guidelines)?
  - [ ] Is rebase and/or squash necessary? If so, please provide the author with appropriate instructions. Also see instructions for [rebase](https://docs.astropy.org/en/latest/development/workflow/development_workflow.html#rebase-if-necessary) and [squash](https://docs.astropy.org/en/latest/development/workflow/development_workflow.html#squash-if-necessary).
  - [ ] Did the CI pass? If no, are the failures related? If you need to run daily and weekly cron jobs as part of the PR, please apply the "Extra CI" label. Codestyle issues can be fixed by the [bot](https://docs.astropy.org/en/latest/development/workflow/development_workflow.html#pre-commit).
  - [ ] Is a change log needed? If yes, did the change log check pass? If no, add the "no-changelog-entry-needed" label. If this is a manual backport, use the "skip-changelog-checks" label unless special changelog handling is necessary.
  - [ ] Is this a big PR that makes a "What's new?" entry worthwhile and if so, is (1) a "what's new" entry included in this PR and (2) the "whatsnew-needed" label applied?
  - [ ] At the time of adding the milestone, if the milestone set requires a backport to release branch(es), apply the appropriate "backport-X.Y.x" label(s) *before* merge.
cmarmo: ... the tests passed locally, I don't understand the reason of the failures  here....
mhvk: The current failures are unrelated - see https://github.com/scientific-python/pytest-doctestplus/pull/242
cmarmo: Sorry for the late answer... life got in the way....

I believe I have added a test checking that the original and copied `_converted` records do not reference the same object...  
cmarmo: @saimn , sorry, I'm lost.... as we discussed offline, I'm trying to access the items of the dictionary `_converted` , but in the particular tested case the dictionary appears to be empty... ? ...
There is at least one scaled column in the table: my understanding was that `_converted` is meant to store the scaled values?
I'm surely missing something... 
pllim: If this is re-milestone to v7.0 or later, please use the new `perf` change log category if applicable, see:

* https://github.com/astropy/astropy/pull/16694
cmarmo: Hello there... I'm sorry I don't know how to check the items of the `_converted` dictionaries as they both appear empty to me.
There is something I don't understand: can someone help? Thanks a lot for your patience.  
saimn: Getting back to this one, the memory leak can be fixed by breaking the references on `_col_weakrefs`:
```
        new = super().copy(order=order)
        new._col_weakrefs = weakref.WeakSet()
        new._coldefs = ColDefs(self._coldefs)
        new._converted = {}
        return new
```
Not sure if resetting `_converted` is needed but I think it's better to be safe. And we would need to test that the BinTableHDU features still work after a copy (e.g. VLA) including for a file read from disk.
cmarmo: @saimn thank you for your patience!
I guess this would go in the 7.x versions, I have updated the changelog entry then.

As you suggested, I have added a test against VLA binary table copies and checked that it was successful on main.
Directly using the [suggested code](https://github.com/astropy/astropy/pull/16143#issuecomment-2442262428) made the test fail with
```
E           OSError: Could not find heap data for the 'QUAL_SPE' variable-length array column.
```

I decided to keep the deepcopy, then, and reset the attributes generating the leak.
Resetting `_converted` is also generating the same `OSError` about heap data, so I only reset `_col_weakrefs`.

saimn: >  Resetting _converted is also generating the same OSError about heap data, so I only reset _col_weakrefs.

Weird but at least we spot that now. I cannot find why so current solution is fine, let's get this fix in.

Patch:
diff --git a/astropy/io/fits/fitsrec.py b/astropy/io/fits/fitsrec.py
index f88a553b00..dedcfafeb1 100644
--- a/astropy/io/fits/fitsrec.py
+++ b/astropy/io/fits/fitsrec.py
@@ -594,6 +594,7 @@ class FITS_rec(np.recarray):
         new = super().copy(order=order)
 
         new.__dict__ = copy.deepcopy(self.__dict__)
+        new._col_weakrefs = weakref.WeakSet()
         return new
 
     @property
diff --git a/docs/changes/io.fits/16143.bugfix.rst b/docs/changes/io.fits/16143.bugfix.rst
new file mode 100644
index 0000000000..8068886d0f
--- /dev/null
+++ b/docs/changes/io.fits/16143.bugfix.rst
@@ -0,0 +1,1 @@
+Fix memory leak in ```BinTableHDU.copy()```


Test patch:
diff --git a/astropy/io/fits/tests/test_table.py b/astropy/io/fits/tests/test_table.py
index a377123843..4457b1cc65 100644
--- a/astropy/io/fits/tests/test_table.py
+++ b/astropy/io/fits/tests/test_table.py
@@ -2858,6 +2858,27 @@ class TestTableFunctions(FitsTestCase):
             t3.teardown_class()
         del t3
 
+    @pytest.mark.skipif(not HAVE_OBJGRAPH, reason="requires objgraph")
+    def test_reference_leak_copyhdu(self):
+        """Regression test for https://github.com/astropy/astropy/issues/15649"""
+
+        def readfile(filename):
+            with fits.open(filename) as hdul:
+                hdu = hdul[1].copy()
+
+        with _refcounting("FITS_rec"):
+            readfile(self.data("memtest.fits"))
+
+    def test_converted_copy(self):
+        """
+        Test that the bintable converted arrays are copied when the fitsrec is.
+        Related to https://github.com/astropy/astropy/issues/15649
+        """
+
+        with fits.open(self.data("memtest.fits")) as hdul:
+            data = hdul[1].data.copy()
+            assert data._converted is not hdul[1].data._converted
+
     def test_dump_overwrite(self):
         with fits.open(self.data("table.fits")) as hdul:
             tbhdu = hdul[1]
@@ -3001,6 +3022,9 @@ class TestVLATables(FitsTestCase):
                 q = toto[1].data.field("QUAL_SPE")
                 assert (q[0][4:8] == np.array([0, 0, 0, 0], dtype=np.uint8)).all()
                 assert toto[1].columns[0].format.endswith("J(1571)")
+                # Test BINTableHDU copy
+                copytoto = toto[1].copy()
+                assert (toto[1].data.field("QUAL_SPE")[0][4:8] == q[0][4:8]).all()
 
         for code in ("PJ()", "QJ()"):
             test(code)