.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "generated/examples/io/skip_create-large-fits.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        Click :ref:`here <sphx_glr_download_generated_examples_io_skip_create-large-fits.py>`
        to download the full example code

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_generated_examples_io_skip_create-large-fits.py:


==========================================
Create a very large FITS file from scratch
==========================================

This example demonstrates how to create a large file (larger than will fit in
memory) from scratch using `astropy.io.fits`.


*By: Erik Bray*

*License: BSD*

.. GENERATED FROM PYTHON SOURCE LINES 17-18

Normally to create a single image FITS file one would do something like:

.. GENERATED FROM PYTHON SOURCE LINES 18-28

.. code-block:: default


    import os

    import numpy as np

    from astropy.io import fits

    data = np.zeros((40000, 40000), dtype=np.float64)
    hdu = fits.PrimaryHDU(data=data)


.. GENERATED FROM PYTHON SOURCE LINES 29-31

Then use the `astropy.io.fits.writeto()` method to write out the new
file to disk

.. GENERATED FROM PYTHON SOURCE LINES 31-34

.. code-block:: default


    hdu.writeto('large.fits')


.. GENERATED FROM PYTHON SOURCE LINES 35-55

However, a 40000 x 40000 array of doubles is nearly twelve gigabytes! Most
systems won't be able to create that in memory just to write out to disk. In
order to create such a large file efficiently requires a little extra work,
and a few assumptions.

First, it is helpful to anticipate about how large (as in, how many keywords)
the header will have in it. FITS headers must be written in 2880 byte
blocks, large enough for 36 keywords per block (including the END keyword in
the final block). Typical headers have somewhere between 1 and 4 blocks,
though sometimes more.

Since the first thing we write to a FITS file is the header, we want to write
enough header blocks so that there is plenty of padding in which to add new
keywords without having to resize the whole file. Say you want the header to
use 4 blocks by default. Then, excluding the END card which Astropy will add
automatically, create the header and pad it out to 36 * 4 cards.

Create a stub array to initialize the HDU; its
exact size is irrelevant, as long as it has the desired number of
dimensions

.. GENERATED FROM PYTHON SOURCE LINES 55-62

.. code-block:: default


    data = np.zeros((100, 100), dtype=np.float64)
    hdu = fits.PrimaryHDU(data=data)
    header = hdu.header
    while len(header) < (36 * 4 - 1):
        header.append()  # Adds a blank card to the end


.. GENERATED FROM PYTHON SOURCE LINES 63-69

Now adjust the NAXISn keywords to the desired size of the array, and write
only the header out to a file. Using the ``hdu.writeto()`` method will cause
astropy to "helpfully" reset the NAXISn keywords to match the size of the
dummy array. That is because it works hard to ensure that only valid FITS
files are written. Instead, we can write just the header to a file using the
`astropy.io.fits.Header.tofile` method:

.. GENERATED FROM PYTHON SOURCE LINES 69-74

.. code-block:: default


    header['NAXIS1'] = 40000
    header['NAXIS2'] = 40000
    header.tofile('large.fits')


.. GENERATED FROM PYTHON SOURCE LINES 75-79

Finally, grow out the end of the file to match the length of the
data (plus the length of the header). This can be done very efficiently on
most systems by seeking past the end of the file and writing a single byte,
like so:

.. GENERATED FROM PYTHON SOURCE LINES 79-90

.. code-block:: default


    with open('large.fits', 'rb+') as fobj:
        # Seek past the length of the header, plus the length of the
        # Data we want to write.
        # 8 is the number of bytes per value, i.e. abs(header['BITPIX'])/8
        # (this example is assuming a 64-bit float)
        # The -1 is to account for the final byte that we are about to
        # write:
        fobj.seek(len(header.tostring()) + (40000 * 40000 * 8) - 1)
        fobj.write(b'\0')


.. GENERATED FROM PYTHON SOURCE LINES 91-92

More generally, this can be written:

.. GENERATED FROM PYTHON SOURCE LINES 92-98

.. code-block:: default


    shape = tuple(header[f'NAXIS{ii}'] for ii in range(1, header['NAXIS']+1))
    with open('large.fits', 'rb+') as fobj:
        fobj.seek(len(header.tostring()) + (np.product(shape) * np.abs(header['BITPIX']//8)) - 1)
        fobj.write(b'\0')


.. GENERATED FROM PYTHON SOURCE LINES 99-108

On modern operating systems this will cause the file (past the header) to be
filled with zeros out to the ~12GB needed to hold a 40000 x 40000 image. On
filesystems that support sparse file creation (most Linux filesystems, but not
the HFS+ filesystem used by most Macs) this is a very fast, efficient
operation. On other systems your mileage may vary.

This isn't the only way to build up a large file, but probably one of the
safest. This method can also be used to create large multi-extension FITS
files, with a little care.

.. GENERATED FROM PYTHON SOURCE LINES 110-111

Finally, we'll remove the file we created:

.. GENERATED FROM PYTHON SOURCE LINES 111-113

.. code-block:: default


    os.remove('large.fits')


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** ( 0 minutes  0.000 seconds)


.. _sphx_glr_download_generated_examples_io_skip_create-large-fits.py:


.. only :: html

 .. container:: sphx-glr-footer
    :class: sphx-glr-footer-example



  .. container:: sphx-glr-download sphx-glr-download-python

     :download:`Download Python source code: skip_create-large-fits.py <skip_create-large-fits.py>`



  .. container:: sphx-glr-download sphx-glr-download-jupyter

     :download:`Download Jupyter notebook: skip_create-large-fits.ipynb <skip_create-large-fits.ipynb>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_