/*****************************************************************************/
// Copyright (c) Illumina 2008
// Author: Richard Shaw
//
// This software is covered by the "Illumina Genome Analyzer Software
// License Agreement" and the "Illumina Source Code License Agreement",
// and certain third party copyright/licenses, and any user of this
// source file is bound by the terms therein (see accompanying files
// Illumina_Genome_Analyzer_Software_License_Agreement.pdf and
// Illumina_Source_Code_License_Agreement.pdf and third party
// copyright/license notices).
/*****************************************************************************/

#include <cassert>
#include <numeric>

#include "Incremental_Accum.h"
#include "Info_Calc.h"

#include "statistics/Info_Contents.h"

/*****************************************************************************/

/*****************************************************************************/

Info_Contents::~Info_Contents()
{
    ;
}

/*****************************************************************************/

void Info_Contents::derive_cumulative_stats()
{
    Incremental_Accum<float> info_inc_accum;
    std::transform(my_info_vec.begin(), my_info_vec.end(),
                   my_cum_info_vec.begin(), info_inc_accum);

    Incremental_Accum<unsigned int> aligned_inc_accum;
    std::transform(my_aligned_vec.begin(), my_aligned_vec.end(),
                   my_cum_aligned_vec.begin(), aligned_inc_accum);

    Incremental_Accum<unsigned int> total_inc_accum;
    std::transform(my_total_vec.begin(), my_total_vec.end(),
                   my_cum_total_vec.begin(), total_inc_accum);
}

/*****************************************************************************/

bool Info_Contents::derive(Base_Calls<unsigned int>& base_calls,
                           std::valarray<unsigned int>& cycle_non_blank_freqs)
{
    typedef Base_Calls<unsigned int> Uint_Base_Calls;
    const Uint_Base_Calls::List base_calls_list
        = base_calls.get_cycle_base_calls_list();
    const unsigned int num_cycles(base_calls_list.size());

    assert(cycle_non_blank_freqs.size() == num_cycles);

    my_info_vec.resize(num_cycles);
    my_aligned_vec.resize(num_cycles);
    my_total_vec.resize(num_cycles);

    my_cum_info_vec.resize(num_cycles);
    my_cum_aligned_vec.resize(num_cycles);
    my_cum_total_vec.resize(num_cycles);

    Cycle_Float_Vec_Iter info_iter(my_info_vec.begin());
    Cycle_Uint_Vec_Iter aligned_iter(my_aligned_vec.begin());
    Cycle_Uint_Vec_Iter total_iter(my_total_vec.begin());
    unsigned int cycle_ind(0);
    const bool exclude_read_blanks(true);
    const bool exclude_ref_blanks(true);

    for (Uint_Base_Calls::List::const_iterator
             base_calls_citer(base_calls_list.begin());
         base_calls_citer != base_calls_list.end();
         ++base_calls_citer, ++info_iter, ++aligned_iter, ++total_iter,
         ++cycle_ind) {
        const Uint_Base_Calls& cycle_base_calls(*base_calls_citer);
        float cycle_info(0.0);
        unsigned int cycle_aligned_count(0);

        Uint_Base_Calls::Const_Range base_calls_range
            = cycle_base_calls.get_crange(Uint_Base_Calls::DIM_REF,
                                          exclude_read_blanks,
                                          exclude_ref_blanks);

        // Iterate over reads.
        for ( ; base_calls_range(); ++base_calls_range) {
            const unsigned int read_base_freq
                = std::accumulate(base_calls_range.begin(),
                                  base_calls_range.end(), 0);
            cycle_aligned_count += read_base_freq;

            if (read_base_freq > 0) {
                Info_Calc info_calc(read_base_freq);


                // Beware the side effects of for_each.
                for (Base_Calls_CIter<unsigned int>
                         freq_citer(base_calls_range.begin());
                     freq_citer != base_calls_range.end(); ++freq_citer) {
                    info_calc(*freq_citer);
                }

                cycle_info += info_calc.info_val();
            }
        }

        *info_iter = cycle_info;
        *aligned_iter = cycle_aligned_count;
        *total_iter = cycle_non_blank_freqs[cycle_ind];
    }

    derive_cumulative_stats();

    return true;
}

/*****************************************************************************/

std::ostream& operator<<(std::ostream& ostrm,
                         const Info_Contents& info_content)
{
    ostrm << "\tBases this cycle\t\tBases so far" << std::endl
          << "Cycle:\tEquiv info:\tAlign:\tTotal:\tEquiv info:\tAlign:\tTotal:"
          << std::endl;

    const size_t num_cycles(info_content.my_info_vec.size());
    const char cycle_prefix_ch(info_content.my_cycle_prefix_ch);

    const unsigned int num_dec_places(2);
    const unsigned int num_pre_pt_places(5);
    const unsigned int width(num_pre_pt_places + 1 + num_dec_places);
    ostrm << std::fixed << std::setprecision(num_dec_places);

    for (size_t cycle_ind(0); cycle_ind < num_cycles; ++cycle_ind) {
        ostrm << cycle_prefix_ch << (cycle_ind + 1) << "\t"
              << std::setw(width)
              << info_content.my_info_vec[cycle_ind] << "\t"
              << info_content.my_aligned_vec[cycle_ind] << "\t"
              << info_content.my_total_vec[cycle_ind] << "\t"
              << std::setw(width)
              << info_content.my_cum_info_vec[cycle_ind] << "\t"
              << info_content.my_cum_aligned_vec[cycle_ind] << "\t"
              << info_content.my_cum_total_vec[cycle_ind] << std::endl;
    }

    return ostrm;
}

/*****************************************************************************/