pca_binning.cc

Go to the documentation of this file.

#include "pca_binning.h"

#include <boost/progress.hpp>

using namespace FullPhysics;
using namespace blitz;

extern "C" {
    void create_bin_uvvswir_v1(int *E_nlayers, int *E_ndat, int *E_maxbins, int *ndat, int *nlay, int *nbin, double *binlims, double *gasdat, double *omega, int *ncnt, int *index);
}

void PCABinning::compute_bins(const SpectralDomain& spec_domain, int channel_index)
{
    Range all = Range::all();

    Array<double, 1> wavenumber(spec_domain.wavenumber());
    boost::shared_ptr<boost::progress_display> progress(new boost::progress_display(wavenumber.rows(), *Logger::stream()));
   
    Array<double, 2> optical_depth(atmosphere->number_layer(), wavenumber.rows(), blitz::ColumnMajorArray<2>());
    Array<double, 2> single_scattering_albedo(optical_depth.shape(), blitz::ColumnMajorArray<2>());

    // Gather optical depth and single scattering albedo for binning routines
    for(int wn_idx = 0; wn_idx < wavenumber.rows(); ++wn_idx) {
        optical_depth(all, wn_idx) = atmosphere->optical_depth_wrt_iv(wavenumber(wn_idx), channel_index).value();
        single_scattering_albedo(all, wn_idx) = atmosphere->single_scattering_albedo_wrt_iv(wavenumber(wn_idx), channel_index).value();

        if (progress) *progress += 1;
    }

    // Call selected binning routine
    int nlayer = atmosphere->number_layer();
    int ndat = wavenumber.rows();
    int nbin = 11;

    binlims_.reference(Array<double, 1>(nbin, blitz::ColumnMajorArray<1>()));
    ncnt_.reference(Array<int, 1>(nbin, blitz::ColumnMajorArray<1>()));
    index_.reference(Array<int, 1>(wavenumber.rows(), blitz::ColumnMajorArray<1>()));
    
    create_bin_uvvswir_v1(&nlayer, &ndat, &nbin, &ndat, &nlayer, &nbin, binlims_.dataFirst(), optical_depth.dataFirst(), single_scattering_albedo.dataFirst(), ncnt_.dataFirst(), index_.dataFirst());

    // Compute binned optical properties, packing into the intermediate variable
}