#ifndef __FILE_LIBMLE_H_SEEN__
#define __FILE_LIBMLE_H_SEEN__ 

/* ----------------------------------------------------------------------------

Copyright (C) 2013.

A. Ronald Gallant

This program is free software; you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation; either version 2 of the License, or
(at your option) any later version.

This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU General Public License for more details.

You should have received a copy of the GNU General Public License along
with this program; if not, write to the Free Software Foundation, Inc.,
51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.

-----------------------------------------------------------------------------*/

#include "libscl.h" 
#include "libmle_base.h"
#include <queue>

namespace libmle {

class mle_mcmc : public mcmc_base {
private:
  proposal_base& T;
  usrmod_base& U;
  INTEGER simulation_size;
  INTEGER stride;
  bool draw_from_posterior;
  REAL temp;
  scl::realmat mode;
  REAL high;
public:
  mle_mcmc(proposal_base& T_fn, usrmod_base& U_mod) 
    : T(T_fn), U(U_mod), simulation_size(1), stride(1),
      draw_from_posterior(true), temp(1.0), mode(), high(-REAL_MAX) { } 
  scl::realmat draw(INT_32BIT& seed, scl::realmat& rho_start,
    scl::realmat& rho_sim, scl::realmat& stats_sim, 
    scl::realmat& pi_sim);
  void set_simulation_size(INTEGER n);
  void set_stride(INTEGER k);
  void set_draw_from_posterior(bool from_posterior);
  void set_temp(REAL temperature);
  void set_mode(const scl::realmat& new_mode, REAL new_high);
  REAL get_temp() { return temp; }
  scl::realmat get_mode() { return mode; }
  REAL get_high() { return high; }
};

class group_move : public proposal_base {
private:
  struct prop_aux {
    REAL    prob;
    scl::realmat Rmat;
    scl::realmat Vinv;
    REAL    scale;
  };
  prop_def pd;
  INTEGER  lrho;
  INTEGER  lpd;
  std::vector<prop_aux> pa;
public:
  group_move(prop_def def);
  scl::den_val operator()(const scl::realmat& rho_old, 
    const scl::realmat& rho_new);
  void draw(INT_32BIT& seed, const scl::realmat& rho_old, 
    scl::realmat& rho_new); 
  INTEGER len_rho() { return lrho; }
  bool transition_is_symmetric() { return true; }
};

class conditional_move : public proposal_base {
private:
  struct possible_move {
    INTEGER y;
    scl::intvec  x;
    REAL    b0;
    scl::realmat b;
    REAL    s;
  };
  INTEGER  lrho;
  INTEGER  lpm;
  std::vector<possible_move> pm;
public:
  conditional_move(prop_def def, std::ostream& detail, bool print);
  scl::den_val operator()(const scl::realmat& rho_old, 
    const scl::realmat& rho_new);
  void draw(INT_32BIT& seed, const scl::realmat& rho_old, 
    scl::realmat& rho_new); 
  INTEGER len_rho() { return lrho; }
};

class minimal_asymptotics : public asymptotics_base {
private:
  const scl::realmat& data;
  INTEGER n_dat;
  mcmc_base& mcmc;
  INTEGER p;
  scl::realmat rho_sum;
  scl::realmat rho_ctr;
  scl::realmat rho_sse;
  scl::realmat mean;
  scl::realmat mode;
  scl::realmat sscp;
  scl::realmat foc;
  INTEGER n_tot;
  REAL high; // temperature does not affect, prior does ifdef ADJUST_FOR_PRIOR
public:
  minimal_asymptotics(const scl::realmat& dat, usrmod_base& mod, 
    mcmc_base& mc) 
    : data(dat), n_dat(dat.get_cols()), mcmc(mc), p(mod.len_rho()), 
      rho_sum(p,1,0.0), rho_ctr(p,1,0.0), rho_sse(p,p,0.0), mean(p,1,0.0), 
      mode(p,1,0.0), sscp(p,p,0.0), foc(p,1,0.0), n_tot(0), high(-REAL_MAX)
    { };
  bool set_asymptotics(const scl::realmat& sim); 
    // Aggregates if called repeatedly.
  void get_asymptotics(scl::realmat& rho_hat, scl::realmat& V_hat, 
    INTEGER& n);
  void get_asymptotics(scl::realmat& rho_mean, scl::realmat& rho_mode, 
    REAL& post_high, scl::realmat& I, scl::realmat& invJ, 
    scl::realmat& foc_hat, INTEGER& reps);
};

class sandwich_asymptotics : public asymptotics_base {
private:
  const scl::realmat& data;
  INTEGER n_dat;
  usrmod_base& model;
  mcmc_base& mcmc;
  INTEGER p;
  scl::realmat rho_sum;
  scl::realmat rho_ctr;
  scl::realmat rho_sse;
  scl::realmat mean;
  scl::realmat mode;
  scl::realmat sscp;
  scl::realmat foc;
  INTEGER n_tot;
  scl::realmat I_mat;
  INTEGER I_reps;
  REAL high; // temperature does not affect, prior does ifdef ADJUST_FOR_PRIOR
  INTEGER lhac;
public:
  sandwich_asymptotics(const scl::realmat& dat, usrmod_base& mod, 
    mcmc_base& mc, INTEGER lags) 
    : data(dat), n_dat(dat.get_cols()), model(mod), mcmc(mc), p(mod.len_rho()), 
      rho_sum(p,1,0.0), rho_ctr(p,1,0.0), rho_sse(p,p,0.0), mean(p,1,0.0), 
      mode(p,1,0.0), sscp(p,p,0.0), foc(p,1,0.0), n_tot(0), I_mat(p,p,0.0), 
      I_reps(0), high(-REAL_MAX), lhac(lags) 
    { };
  bool set_asymptotics(const scl::realmat& sim); 
    // Aggregates if called repeatedly.
  void get_asymptotics(scl::realmat& rho_hat, scl::realmat& V_hat, 
    INTEGER& n);
  void get_asymptotics(scl::realmat& rho_mean, scl::realmat& rho_mode, 
    REAL& post_high, scl::realmat& I, scl::realmat& invJ, 
    scl::realmat& foc_hat, INTEGER& reps);
};

}

#endif
