merge.cpp

Go to the documentation of this file.

/**
 * @file
 * @ingroup daq_dpm_merge
 * @copyright ESO - European Southern Observatory
 */
#include "merge.hpp"
 
#include <fmt/format.h>
#include <log4cplus/loggingmacros.h>
 
#include <daq/fits/cfitsio.hpp>
 
namespace daq::dpm::merge {
 
void CheckUnmergedPrimaryHduData(log4cplus::Logger const& logger,
                                 StatusReporter& status_reporter,
                                 FitsFileSource const& source) {
    int status = 0;
    int type = 0;
    fitsfile* file = source.GetFitsFile();
 
    // Select primary HDU (1)
    (void)fits_movabs_hdu(file, 1, &type, &status);
    // As primary HDU must be empty or contain image (and not tables) we just check the image axis
    // dimensions.
    auto const max_dim = 10;
    long naxes[max_dim] = {0};
    int naxis = 0;
    fits_get_img_dim(file, &naxis, &status);
    fits_get_img_size(file, max_dim, naxes, &status);
    if (status != 0) {
        // Non-conforming FITS -> ignore.
        return;
    }
 
    naxis = naxis < max_dim ? naxis : max_dim;
    bool has_data = false;
    for (auto axis = 0; axis < naxis; ++axis) {
        if (naxes[axis] > 0u) {
            // non-empty axis
            has_data = true;
            break;
        }
    }
 
    if (has_data) {
        auto dim_str = std::string();
        for (auto axis = 0; axis < naxis; ++axis) {
            if (naxes[axis] > 0u) {
                auto dim = naxes[axis];
                if (axis > 0) {
                    dim_str += fmt::format(", {}", dim);
                } else {
                    dim_str += fmt::format("{}", dim);
                }
            }
        }
        auto msg = fmt::format(
            "Unmerged data: Source file has data in primary HDU that cannot be merged: '{}' "
            "produced file '{}' that was copied from '{}' has dimensions ({})",
            source.GetName(),
            source.GetFilePath().c_str(),
            source.GetLocation(),
            dim_str);
        LOG4CPLUS_WARN(logger, msg);
        if (source.AlertUnmergeable()) {
            status_reporter.PostAlert(
                fmt::format(
                    "unmerged-{}-{}", source.GetName(), source.GetFilePath().filename().native()),
                msg);
        }
    }
}
 
fits::KeywordVector CompileKeywords(Operations ops, FitsKeywordsSource const& source) {
    try {
        return source.GetKeywordRuleProcessor().Process(
            source.GetKeywords(),
            source.GetInitialKeywords().value_or(KeywordRuleProcessor::DefaultRule::User));
    } catch (...) {
        auto msg = fmt::format("Processing keyword rules failed for JSON keywords source '{}'",
                               source.GetName());
        LOG4CPLUS_ERROR(ops.logger, msg);
        std::throw_with_nested(std::runtime_error(msg));
    }
}
 
fits::KeywordVector CompileKeywords(Operations ops, FitsFileSource const& source, int hdu_num) {
    auto logger = log4cplus::Logger::getInstance("daq.dpmmerge");
    auto literal_kws = fits::ReadKeywords(source.GetFitsFile(), hdu_num);
    fits::KeywordVector result;
    result.reserve(literal_kws.size());
    std::copy(std::make_move_iterator(literal_kws.begin()),
              std::make_move_iterator(literal_kws.end()),
              std::back_inserter(result));
 
    try {
        return source.GetKeywordRuleProcessor().Process(
            result, source.GetInitialKeywords().value_or(KeywordRuleProcessor::DefaultRule::User));
    } catch (...) {
        auto msg = fmt::format("Processing keyword rules failed for FITS file keywords source '{}'",
                               source.GetFilePath().string());
        LOG4CPLUS_ERROR(logger, msg);
        std::throw_with_nested(std::runtime_error(msg));
    }
}
 
/**
 * Compiles unique keywords and formats primary HDU keywords
 */
fits::KeywordVector CompilePrimaryHduKeywords(Operations ops,
                                              TargetSource& target,
                                              std::vector<SourceTypes> const& sources) {
    auto logger = log4cplus::Logger::getInstance("daq.dpmmerge");
    constexpr const int primary_hdu_num = 1;
    fits::KeywordVector result;
    auto literal_kws = fits::ReadKeywords(target.GetFitsFile(), primary_hdu_num);
    LOG4CPLUS_DEBUG(logger, "Read keywords from " << target.GetFilePath());
    std::for_each(std::begin(literal_kws), std::end(literal_kws), [&](auto const& kw) {
        LOG4CPLUS_DEBUG(logger, kw);
    });
    result.reserve(literal_kws.size());
    std::copy(std::begin(literal_kws), std::end(literal_kws), std::back_inserter(result));
    // Filter
    try {
        // note: By default we keep all keywords in the target
        result = target.GetKeywordRuleProcessor().Process(
            result, target.GetInitialKeywords().value_or(KeywordRuleProcessor::DefaultRule::All));
        LOG4CPLUS_DEBUG(logger, "Result after keyword processing: " << target.GetFilePath());
        std::for_each(std::begin(result), std::end(result), [&](auto const& kw) {
            LOG4CPLUS_DEBUG(logger, kw);
        });
    } catch (...) {
        auto msg = fmt::format("Processing keyword rules failed for target FITS file '{}'",
                               target.GetFilePath().string());
        LOG4CPLUS_ERROR(logger, msg);
        std::throw_with_nested(std::runtime_error(msg));
    }
 
    for (auto const& source : sources) {
        std::string path_or_name;
        auto kws = std::visit(
            [&](auto const& source) -> fits::KeywordVector {
                using T = std::decay_t<decltype(source)>;
                if constexpr (std::is_same_v<T, FitsKeywordsSource>) {
                    path_or_name = fmt::format("{}: (keyword list)", source.GetName());
                    return CompileKeywords(ops, source);
                } else if constexpr (std::is_same_v<T, FitsFileSource>) {
                    path_or_name =
                        fmt::format("{}: {}", source.GetName(), source.GetFilePath().native());
                    return CompileKeywords(ops, source, 1);
                }
            },
            source);
        LOG4CPLUS_DEBUG(logger, "Updating with keywords from : " << path_or_name);
        std::for_each(
            std::begin(kws), std::end(kws), [&](auto const& kw) { LOG4CPLUS_DEBUG(logger, kw); });
        // Don't overwrite conflicting keywords
        fits::UpdateKeywords(result, kws, fits::ConflictPolicy::Skip);
    }
 
    return result;
}
 
std::vector<fits::LiteralKeyword> FormatKeywords(fits::KeywordVector::const_iterator begin,
                                                 fits::KeywordVector::const_iterator end,
                                                 KeywordFormatter& fmt) {
    auto const& logger = log4cplus::Logger::getInstance("daq.dpmmerge");
    std::vector<fits::LiteralKeyword> result;
    LOG4CPLUS_DEBUG(logger, "Formatting keywords ...");
    std::transform(
        begin, end, std::back_inserter(result), [&](auto const& kw) -> fits::LiteralKeyword {
            LOG4CPLUS_DEBUG(logger, "Formatting keyword input:  \"" << kw << "\"");
            auto formatted = fmt.FormatKeyword(kw);
            LOG4CPLUS_DEBUG(logger, "Formatting keyword result: \"" << formatted << "\"");
            return formatted;
        });
    LOG4CPLUS_DEBUG(logger, "Formatting keywords done.");
    return result;
}
 
template <class Container>
void LogKeywords(log4cplus::Logger const& logger, Container const& keywords) {
    std::for_each(std::begin(keywords), std::end(keywords), [&](auto const& kw) {
        LOG4CPLUS_DEBUG(logger, kw);
    });
}
 
void MergePrimaryHduKeywords(Operations ops,
                             Params const& params,
                             TargetSource& target,
                             std::vector<SourceTypes> const& sources,
                             bool dry_run) {
    auto const& logger = ops.logger;
    LOG4CPLUS_INFO(logger, "Merge primary HDU keywords");
    LOG4CPLUS_INFO(logger, "Compile primary keywords");
    constexpr const int primary_hdu_num = 1;
    auto primary_hdu_keywords = CompilePrimaryHduKeywords(ops, target, sources);
 
    // Add ORIGFILE and ARCFILE last (will be sorted to be last of value keywords)
    fits::KeywordVector mandatory;
    mandatory.emplace_back(std::in_place_type<fits::ValueKeyword>, "ORIGFILE", params.origfile);
    mandatory.emplace_back(std::in_place_type<fits::ValueKeyword>, "ARCFILE", params.arcfile);
    fits::InsertKeywords(
        primary_hdu_keywords, primary_hdu_keywords.end(), mandatory.begin(), mandatory.end());
 
    LOG4CPLUS_INFO(logger, "Format keywords");
    auto formatted = FormatKeywords(
        std::begin(primary_hdu_keywords), std::end(primary_hdu_keywords), ops.keyword_formatter);
    {
        LOG4CPLUS_INFO(logger, "Sort keywords");
        ops.keyword_sorter.SortKeywords(formatted);
        LOG4CPLUS_DEBUG(logger, "Sorted keywords");
        LogKeywords(logger, formatted);
    }
 
    // Write keywords back
    // Write keywords will make room as necessary.
    if (!dry_run) {
        LOG4CPLUS_INFO(logger, "Clear keywords to make room for writing back sorted keywords.");
        fits::DeleteAllKeywords(target.GetFitsFile(), primary_hdu_num);
        LOG4CPLUS_INFO(logger, "Writing keywords");
        std::optional<ssize_t> remaining_size;
        fits::WriteKeywords(target.GetFitsFile(), primary_hdu_num, formatted, &remaining_size);
        if (remaining_size) {
            if (*remaining_size < 0) {
                auto needed = -*remaining_size;
                // Reallocation occurred
                auto msg = fmt::format(
                    "Writing keywords required resizing of primary HDU: Add space for at least "
                    "{} "
                    "keywords to avoid resize",
                    needed);
                ops.status_reporter.PostAlert("primary_hdu_resize", msg);
            } else {
                auto msg = fmt::format("Primary HDU keyword space remaining: {} ", *remaining_size);
                LOG4CPLUS_INFO(logger, msg);
            }
        }
    } else {
        LOG4CPLUS_INFO(logger, "Writing keywords SKIPPED (dry-run)");
    }
}
 
/**
 * Copy all extensions from source to target
 */
void CopyExtensions(Operations ops,
                    TargetSource& target,
                    FitsFileSource const& source,
                    bool dry_run) {
    auto const& logger = ops.logger;
    LOG4CPLUS_INFO(
        logger,
        "Merging HDU extensions from " << source.GetName() << "(" << source.GetFilePath() << ")");
 
    // Check for unmerged primary HDU data.
    CheckUnmergedPrimaryHduData(ops.logger, ops.status_reporter, source);
 
    int status = 0;
    fitsfile* source_fits = source.GetFitsFile();
    fitsfile* target_fits = target.GetFitsFile();
 
    int num_hdus = 0;
    fits_get_num_hdus(source_fits, &num_hdus, &status);
    if (status != 0) {
        auto msg =
            fmt::format("Failed to get number of HDUs from '{}'", target.GetFilePath().c_str());
        LOG4CPLUS_ERROR(logger, msg);
        throw fits::CfitsioError(status, msg);
    }
 
    if (num_hdus == 1) {
        // @todo There should be a proper check of the first HDU to see that it does not contain
        // any data!
        LOG4CPLUS_INFO(ops.logger,
                       "Note: No HDU extensions to merge from " << source.GetFilePath());
        return;
    }
 
    // Select HDU 2 to copy
    fits::SelectHduNum(source_fits, 2);
    if (!dry_run) {
        // @todo Execute HDU by HDU to allow aborting in-between?
        int previous = false;
        int current = true;
        int following = true;
        fits_copy_file(source_fits, target_fits, previous, current, following, &status);
        if (status != 0) {
            auto const msg = "FITS function fits_copy_file failed";
            LOG4CPLUS_ERROR(logger, msg);
            throw fits::CfitsioError(status, msg);
        }
    } else {
        LOG4CPLUS_INFO(logger,
                       "Merging HDU extensions from " << source.GetName() << "("
                                                      << source.GetFilePath()
                                                      << ") SKIPPED (dry-run)");
    }
}
 
void MergeHduExtensions(Operations ops,
                        TargetSource& target,
                        std::vector<SourceTypes> const& sources,
                        bool dry_run) {
    auto const& logger = ops.logger;
    LOG4CPLUS_INFO(logger, "Merging HDU extensions");
 
    for (auto const& source_var : sources) {
        if (!std::holds_alternative<FitsFileSource>(source_var)) {
            continue;
        }
        FitsFileSource const& source = std::get<FitsFileSource>(source_var);
        try {
            CopyExtensions(ops, target, source, dry_run);
        } catch (...) {
            std::throw_with_nested(
                std::runtime_error(fmt::format("Failed to copy HDU extensions from '{}' to '{}'",
                                               source.GetFilePath().c_str(),
                                               target.GetFilePath().c_str())));
        }
    }
    LOG4CPLUS_INFO(logger, "Merging HDU extensions completed successfully");
}
 
void UpdateChecksums(Operations ops, TargetSource& target, bool dry_run) {
    auto const& logger = ops.logger;
    LOG4CPLUS_INFO(logger, "Updating checksums for all HDUs");
    int status = 0;
    fitsfile* target_fits = target.GetFitsFile();
    int num_hdus = 0;
    fits_get_num_hdus(target_fits, &num_hdus, &status);
    if (status != 0) {
        auto msg =
            fmt::format("Failed to get number of HDUs from '{}'", target.GetFilePath().c_str());
        LOG4CPLUS_ERROR(logger, msg);
        throw fits::CfitsioError(status, msg);
    }
 
    if (!dry_run) {
        for (int hdu_num = 1; hdu_num <= num_hdus; ++hdu_num) {
            LOG4CPLUS_DEBUG(logger, "Updating checksum for HDU " << hdu_num);
            fits::WriteChecksum(target_fits, hdu_num);
        }
    } else {
        LOG4CPLUS_INFO(logger, "Updating checkum keywords SKIPPED (dry-run)");
    }
    LOG4CPLUS_INFO(logger, "Updating checkum keywords completed for all HDUs successfully");
}
 
void Merge(Operations ops,
           Params const& params,
           TargetSource& target,
           std::vector<SourceTypes> const& sources,
           bool dry_run) {
    auto const& logger = ops.logger;
    /*
     * 1. Merge target primary HDU keywords
     *    Since adding keywords to primary HDU may shift the data section we do this first,
     * before copying any HDU extensions
     *
     *    1. Compile keywords to be merged to target primary HDU.
     *    2. Format and validate keywords.
     *       - Already formatted keywords (literal keywords) are compared with expected format.
     *       - Unformatted keywords are formatted using dictionary format string.
     *    3. Sort keywords
     *    4. Delete all existing keywords.
     *    4. Write keywords cleanly (remove all then write back)
     * 2. Copy HDU extensions to target.
     */
    // @todo Validate that there is no primary HDU data in any of the sources before starting
    // the process.
    LOG4CPLUS_INFO(logger, "Starting merge operation");
    try {
        MergePrimaryHduKeywords(ops, params, target, sources, dry_run);
    } catch (...) {
        std::throw_with_nested(std::runtime_error("Failed to merge primary HDU keywords"));
    }
 
    try {
        MergeHduExtensions(ops, target, sources, dry_run);
    } catch (...) {
        std::throw_with_nested(std::runtime_error("Failed to merge primary HDU keywords"));
    }
    try {
        UpdateChecksums(ops, target, dry_run);
    } catch (...) {
        std::throw_with_nested(std::runtime_error("Failed to update checksums"));
    }
 
    LOG4CPLUS_INFO(logger, "Completed successfully");
}
 
}  // namespace daq::dpm::merge