When dealing with large datasets (images or 3D datacubes), it is very memory and time consuming to read the whole dataset into memory before operations. Nearly every time consuming operation in Gnuastro (not ConvertType for example) can be done on a small part of the input data independently. Fortunately CFITSIO can also read and write from and to parts of a file. So I would like to suggest that we design a highly scalabe input/output+operation system, where for example huge datasets (much larger than the RAM) can be easily processed: for example in warping or to run NoiseChisel on.

Fortunately a very rudimentary first step is already in place for the operations of NoiseChisel and SubtractSky: mesh.[ch]. But all these functions are built for use within RAM. Here I am suggesting to extend the definition and implementation in such a way that it merges with operations in fitsarrayvv.[ch] so when a program thread (CPU or GPU, see task #13787) needs only one part of an image, it only reads that one part, does its operations and then saves its results to the proper part of an output file. This way the RAM will never get full and parallel (CPU or GPU) programming can create a huge speed boost.

The proposed mesh structure should be able to allow a certain amount of margin (or extra pixels) for each mesh, so operations like spatial convolution or warping can also be done in this embarrassingly parallel manner and without relying on a huge RAM. Particularly with SSD devices that are becomming popular, the read/write speed is not a major issue any more.

Like NoiseChisel and SubtractSky, the user can set the mesh parameters, so if they don't have any RAM issue, they can simply only define one mesh and they won't have to worry about read/write speeds slowing them down. However, when the datasets are huge (like a MUSE datacube, or one HSC or LSST exposure), this is very important.

If it becomes too complex, this mesh structure can be made over an image in a separate utility and all the mesh parameters can be added as headers in the FITS file. Then any program that needs them will read those headers and use them as a basis to operate on the image in a highly scalable way.

An array can also be defined to keep one value per mesh in operations like finding the Sky value and its standard deviation in NoiseChisel, SubtractSky and MakeCatalog.

These are all just raw suggestions and are a more generalized case for task #13759 (reducing the memory usage of MakeCatalog). Certainly as the implementaion starts, the details can/will change from what is suggested here. The main goal however is to allow a high level of scalability (more parallel CPU and GPU operations and reducing the dependancy on the RAM)