secorder_input.hpp 2.27 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93
#ifndef SOINPUT_HPP
#define SOINPUT_HPP
#include<iostream>
#include<gsl/gsl_math.h>
#include<gsl/gsl_cdf.h>
#include<gsl/gsl_vector.h>
#include<gsl/gsl_matrix.h>
#include<gsl/gsl_rng.h>
#include<gsl/gsl_randist.h>
#include "calc_rhos.hpp"

class secorder_input  {
  int N_nodes[2];
  double input_p0, input_p1;
  double input_corr0, input_corr1, input_cc_01, ;
  double rho00, rho01, rho11;
  int status;
  double sigma0, sigma1;
  
  double sqrt_diag[2];
  double sqrt_offdiag[2][2];

  gsl_matrix *input_trains;
  double dt_batch, dt_internal;
  int N_intervals_batch;
  double input_rate0, input_rate1;

  int setup_params();

public:
  // secorder_input(double input_p_, double input_corr_, int N_nodes_) {
  //   set_stats(input_p_, input_corr_, N_nodes_);
  // }
  secorder_input() {
    input_trains=0;
  }

  secorder_input(double input_rate0_, double input_rate1_, 
		 double input_corr0_, double input_corr1_,
		 double input_cc_01_, double dt_, 
		 int N_nodes0_, int N_nodes1_) {
    input_trains=0;
    setup_trains(input_rate0_, input_rate1_,
		 input_corr0_,  input_corr1_,  input_cc_01_, 
		 dt_, N_nodes0_, N_nodes1_);
  }
  ~secorder_input() {
    if(input_trains)
      gsl_matrix_free(input_trains);
  }

  int set_stats(double input_p0_, double input_p1_, 
		double input_corr0_, double input_corr1_,
		double input_cc_01_,
		int N_nodes0_=0, int N_nodes1_=0) {
    if (input_p0_>1) {
	std::cerr<<"secorder_input::set_stats: input probability 0 > 1 !!"<<std::endl;
	exit(-1);
    }
    if (input_p1_>1) {
	std::cerr<<"secorder_input::set_stats: input probability 1 > 1 !!"<<std::endl;
	exit(-1);
    }
    input_p0 = input_p0_;
    input_p1 = input_p1_;
    input_corr0 = input_corr0_;
    input_corr1 = input_corr1_;
    input_cc_01 = input_cc_01_;
    if(N_nodes0_) 
      N_nodes[0] = N_nodes0_;
    if(N_nodes1_) 
      N_nodes[1] = N_nodes1_;
    status=setup_params();

    if(status) {
      std::cerr << "secorder_input::set_stats: failed" << std::endl;
      exit(-1);
    }
    return 0;
  }

  int sample(gsl_vector *inputs, gsl_rng *rng);

  int setup_trains(double input_rate0_, double input_rate1_, 
		   double input_corr0_, double input_corr1_,
		   double input_cc_01_, double dt_, 
		   int N_nodes0_, int N_nodes1_);

  gsl_matrix *sample_trains(double time, gsl_rng *rng);
};


#endif