<commit_before>#ifndef __FAST5_SUMMARY_HPP
#define __FAST5_SUMMARY_HPP

#include <array>
#include <string>
#include <vector>

#ifndef H5_HAVE_THREADSAFE
#include <mutex>
#endif

#include "Pore_Model.hpp"
#include "fast5.hpp"
#include "alg.hpp"

template < typename Float_Type = float >
class Fast5_Summary
{
public:
    typedef Pore_Model< Float_Type > Pore_Model_Type;
    typedef Pore_Model_Parameters< Float_Type > Pore_Model_Parameters_Type;
    typedef Event< Float_Type > Event_Type;
    typedef Event_Sequence< Float_Type > Event_Sequence_Type;

    std::string file_name;
    std::string base_file_name;
    std::string read_id;
    std::array< std::string, 2 > preferred_model;
    std::array< std::map< std::string, Pore_Model_Parameters_Type >, 2 > params;
    std::array< unsigned, 4 > strand_bounds;
    unsigned num_ed_events;
    float abasic_level;
    bool have_ed_events;
    bool valid;

    static unsigned& min_read_len()
    {
        static unsigned _min_read_len = 1000;
        return _min_read_len;
    }

    // from fast5 file
    std::vector< fast5::EventDetection_Event_Entry > ed_events;
    // filtered
    std::array< Event_Sequence_Type, 2 > events;

    Fast5_Summary() : valid(false) {}
    Fast5_Summary(const std::string fn, const Pore_Model_Dict_Type& models) : valid(false) { summarize(fn, models); }

    void summarize(const std::string& fn, const Pore_Model_Dict_Type& models)
    {
        valid = true;
        file_name = fn;
        auto pos = file_name.find_last_of('/');
        base_file_name = (pos != std::string::npos? file_name.substr(pos + 1) : file_name);
        if (base_file_name.substr(base_file_name.size() - 6) == ".fast5")
        {
            base_file_name.resize(base_file_name.size() - 6);
        }
        fast5::File f(file_name);
        have_ed_events = f.have_eventdetection_events();
        num_ed_events = 0;
        abasic_level = 0.0;
        if (have_ed_events)
        {
            load_ed_events(f);
            auto ed_params = f.get_eventdetection_event_parameters();
            num_ed_events = ed_events.size();
            read_id = ed_params.read_id;
            if (read_id.empty())
            {
                read_id = base_file_name;
            }
            abasic_level = detect_abasic_level();
            if (abasic_level > 1.0)
            {
                detect_strands();
            }
            // compute initial model scalings
            load_events();
            for (unsigned st = 0; st < 2; ++st)
            {
                if (events[st].size() < min_read_len())
                {
                    continue;
                }
                auto r = alg::mean_stdv_of< Float_Type >(events[st], [] (const Event_Type& ev) { return ev.mean; });
                for (const auto& p : models)
                {
                    if (p.second.strand() == st or p.second.strand() == 2)
                    {
                        Pore_Model_Parameters_Type param;
                        param.scale = r.second / p.second.stdv();
                        param.shift = r.first - param.scale * p.second.mean();
                        LOG("Fast5_Summary", debug)
                            << "read [" << read_id << "] strand [" << st
                            << "] model [" << p.first << "] initial parameters "
                            << param << std::endl;
                        params[st][p.first] = std::move(param);
                    }
                }
            }
            ed_events.clear();
        }
    }

    void load_events()
    {
        assert(valid and have_ed_events);
        if (ed_events.empty())
        {
#ifndef H5_HAVE_THREADSAFE
            static std::mutex fast5_mutex;
            std::lock_guard< std::mutex > fast5_lock(fast5_mutex);
#endif
            fast5::File f(file_name);
            load_ed_events(f);
            f.close();<commit_msg>fix: crash if number of ed events too small<commit_after>#ifndef __FAST5_SUMMARY_HPP
#define __FAST5_SUMMARY_HPP

#include <array>
#include <string>
#include <vector>

#ifndef H5_HAVE_THREADSAFE
#include <mutex>
#endif

#include "Pore_Model.hpp"
#include "fast5.hpp"
#include "alg.hpp"

template < typename Float_Type = float >
class Fast5_Summary
{
public:
    typedef Pore_Model< Float_Type > Pore_Model_Type;
    typedef Pore_Model_Parameters< Float_Type > Pore_Model_Parameters_Type;
    typedef Event< Float_Type > Event_Type;
    typedef Event_Sequence< Float_Type > Event_Sequence_Type;

    std::string file_name;
    std::string base_file_name;
    std::string read_id;
    std::array< std::string, 2 > preferred_model;
    std::array< std::map< std::string, Pore_Model_Parameters_Type >, 2 > params;
    std::array< unsigned, 4 > strand_bounds;
    unsigned num_ed_events;
    float abasic_level;
    bool have_ed_events;
    bool valid;

    static unsigned& min_read_len()
    {
        static unsigned _min_read_len = 1000;
        return _min_read_len;
    }

    // from fast5 file
    std::vector< fast5::EventDetection_Event_Entry > ed_events;
    // filtered
    std::array< Event_Sequence_Type, 2 > events;

    Fast5_Summary() : valid(false) {}
    Fast5_Summary(const std::string fn, const Pore_Model_Dict_Type& models) : valid(false) { summarize(fn, models); }

    void summarize(const std::string& fn, const Pore_Model_Dict_Type& models)
    {
        valid = true;
        file_name = fn;
        auto pos = file_name.find_last_of('/');
        base_file_name = (pos != std::string::npos? file_name.substr(pos + 1) : file_name);
        strand_bounds = { 0, 0, 0, 0 };
        if (base_file_name.substr(base_file_name.size() - 6) == ".fast5")
        {
            base_file_name.resize(base_file_name.size() - 6);
        }
        fast5::File f(file_name);
        have_ed_events = f.have_eventdetection_events();
        num_ed_events = 0;
        abasic_level = 0.0;
        if (have_ed_events)
        {
            load_ed_events(f);
            auto ed_params = f.get_eventdetection_event_parameters();
            num_ed_events = ed_events.size();
            read_id = ed_params.read_id;
            if (read_id.empty())
            {
                read_id = base_file_name;
            }
            abasic_level = detect_abasic_level();
            if (abasic_level > 1.0)
            {
                // compute initial model scalings
                detect_strands();
                load_events();
                for (unsigned st = 0; st < 2; ++st)
                {
                    if (events[st].size() < min_read_len())
                    {
                        continue;
                    }
                    auto r = alg::mean_stdv_of< Float_Type >(events[st], [] (const Event_Type& ev) { return ev.mean; });
                    for (const auto& p : models)
                    {
                        if (p.second.strand() == st or p.second.strand() == 2)
                        {
                            Pore_Model_Parameters_Type param;
                            param.scale = r.second / p.second.stdv();
                            param.shift = r.first - param.scale * p.second.mean();
                            LOG("Fast5_Summary", debug)
                                << "read [" << read_id << "] strand [" << st
                                << "] model [" << p.first << "] initial parameters "
                                << param << std::endl;
                            params[st][p.first] = std::move(param);
                        }
                    }
                }
                ed_events.clear();
            } // if abasic_level > 1.0
        } // if have_ed_events
    }

    void load_events()
    {
        assert(valid and have_ed_events);
        if (ed_events.empty())
        {
#ifndef H5_HAVE_THREADSAFE
            static std::mutex fast5_mutex;
            std::lock_guard< std::mutex > fast5_lock(fast5_mutex);
#endif
            fast5::File f(file_name);
            load_ed_events(f);
            f.close();
{'fraction_non_alphanumeric': 0.11075161772025884, 'fraction_numerical': 0.009706321553011448, 'mean_word_length': 1.812106368089573, 'pattern_counts': {'":': 0, '<': 77, '<?xml version=': 0, '>': 40, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '643292', 'n_tokens_mistral': 2538, 'n_tokens_neox': 2400, 'n_words': 493}