
.. 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 18-19

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

.. GENERATED FROM PYTHON SOURCE LINES 19-26

.. 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 27-29

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

.. GENERATED FROM PYTHON SOURCE LINES 29-32

.. code-block:: default


    hdu.writeto('large.fits')


.. GENERATED FROM PYTHON SOURCE LINES 33-53

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 53-60

.. 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 61-67

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 67-72

.. code-block:: default


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


.. GENERATED FROM PYTHON SOURCE LINES 73-77

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 77-88

.. 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 89-90

More generally, this can be written:

.. GENERATED FROM PYTHON SOURCE LINES 90-96

.. code-block:: default


    shape = tuple(header['NAXIS{0}'.format(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 97-106

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 108-109

Finally, we'll remove the file we created:

.. GENERATED FROM PYTHON SOURCE LINES 109-111

.. 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>`_
