#include <sli/fitscc.h>
#include <sli/mdarray_statistics.h>
#include <stdlib.h>

using namespace sli;


int main(int argc, char *argv[])
{
    if (argc < 3) {
        sli__eprintf("usage: %s OUTPUT FLAT1 FLAT2 ...\n", argv[0]);
        exit(1);
    }

    // コマンド引数とargv, argcの関係
    // ./combine-float out.fits flat1.fits flat2.fits flat3.fits
    // |               |        |          |          |
    // argv[0]         argv[1]  argv[2]    argv[3]    argv[4]
    // argc == 5

    const int n_flat = argc - 2;

    fitscc fits;
    mdarray_float stack;

    for (int i = 0;  i < n_flat;  i++) {
        const char *flat_file = argv[i + 2];
        fits.read_stream(flat_file);

        const mdarray_float &data = fits.image(0L).float_array();

        if (i == 0)
            stack.resize_3d(data.length(0), data.length(1), n_flat);

        stack.paste(data / md_median(data), 0, 0, i);
    }

    stack.dprint(); // stackの内容を確認
    fits.image(0L).float_array() = md_median_small_z(stack);

    fits.write_stream(argv[1]);

    return 0;
}