#include <stdlib.h>
#include "fileutils.h"
#include "math.h"
#include "unistd.h"

#include "check_stdp_freq-dep.h"
#include "fileutils.cpp"
#include "model_switch.h"

using namespace std;

int main(int argc, char **argv) {
  // check cmd line sanity
  if (argc != 5) {
    fprintf(stderr, "Wrong argument count\n\n"
	    "Call format:\n"
	    "%s\n\t"
	    "performance out\n\t"
	    "trace cmd out\n\t"
	    "global out\n\t"
	    "spike out\n\t"
	    "\n"
	    "Special names allowed:\n\t- (standart input)\n\t0 (/dev/null)\n", argv[0]);
    return -1;
  }

  Trainer *t = new Trainer(argc, argv);
  t->run();

  pthread_join(t->thread_write, NULL);
}

Trainer::Trainer(int argc, char** argv) {
  // init vars
  currentEpoch  = 0;
  epochDuration  = 10.0; // [s]
  neurons        = 1000; // number of neurons to send noise to
  voltage        = 0.1;  // [V]
  md             = 1.2; // 10.0;
  mss            = 1.1; //1.2;
  fs             = 100;   // number of frequencies to try
  frd            = 1.0;  // relative difference between to frequencies (f_i+1 = frd * f_i)
  fad            = 0.5;  // absolute difference between to frequencies (f_i+1 = fad + f_i)
  
  // open all file descriptors in an order complementary to the simulators one
  // to avoid deadlocks
  fd_spike_out       = fd_magic(argv[4], true);
  fd_global_out      = fd_magic(argv[3], true);
  fd_performance_out = fd_magic(argv[1], true);
  fd_trace_out       = fd_magic(argv[2], true);

  // create read and write threads
  pthread_create(&thread_write, NULL, (void* (*)(void*)) &write_spikes, this);
}

void Trainer::run() {
  char *str_trace = "%f; synapse\n";

  // init global sim variables
  MS_Global msg;
  msg_init(msg);
  msg.dopamin_level = 0.0;
  
  double ta = 0.009821, //0.0088541,
         la = 0.140249; // 0.126445;

  /*
  // loop over both vars to examine
  for (msg.stdp_tau_plus  = msg.stdp_tau_minus / md;
       msg.stdp_tau_plus <= msg.stdp_tau_minus * md;
       msg.stdp_tau_plus *= mss) {
    for (msg.stdp_lambda_plus  = msg.stdp_lambda_minus / md;
	 msg.stdp_lambda_plus <= msg.stdp_lambda_minus * md;
	 msg.stdp_lambda_plus *= mss) {*/
  // loop over both vars to examine
  for (msg.stdp_tau_plus  = ta / md;
       msg.stdp_tau_plus <= ta * md;
       msg.stdp_tau_plus *= mss) {
    for (msg.stdp_lambda_plus  = la / md;
	 msg.stdp_lambda_plus <= la * md;
	 msg.stdp_lambda_plus *= mss) {
      
      // print the parameters to the performance output
      msg_print(msg, fd_performance_out);
      fprintf(fd_performance_out, "\n");

      // print the global params
      fprintf(fd_global_out, "%f, ", currentEpoch * epochDuration);
      msg_print(msg, fd_global_out);
      fprintf(fd_global_out, "\n");

      // let the simulation proceed
      fprintf(fd_trace_out, str_trace, epochDuration);
      currentEpoch++;
 
      // repeat this 2*n-1 times (n=number of different frequency trials)
      for (int i=0; i < 2*fs-1; i++) {
	fprintf(fd_trace_out, "\n");
	currentEpoch++;
      }
    }
  }

  fclose(fd_trace_out);
  fclose(fd_global_out);
}

// ---- send indepenent poisson noise w/ increasing fequency----
void *write_spikes(Trainer *t) {
  // calculate how often we have to try all frequencies (=outer loop)
  // WARN: ignore minor numerical instabilities
  int max = (int) floor(2.0 * log(t->md) / log(t->mss) ) + 1;
  max *= max; // there are two nested loops of the same size

  double time = 0.0; // global time (that one send to the simulator)

  // for each paramter config (set in the main routine)
  for (int i=0; i<max; i++) {

    double freq = 1.0;

    // examine a set of frequencies
    for (int j=0; j < t->fs; j++) {
      // send out the spikes
      double localtime = 0.0;
      double nextRefSpike = 0.0;
      double refFreq = 10.0; // [Hz]
      int dst = -1;
      while (localtime < t->epochDuration) {
	// starting with the second call ...
	if (dst != -1) {
	  // check if we have to send a spike to the ref neuron
	  if (localtime > nextRefSpike) {
	    fprintf(t->fd_spike_out, "%f, %d, %f\n", time + nextRefSpike, 0, t->voltage);
	    nextRefSpike += 1.0 / refFreq;
	  }

	  // send spike to the simulator
	  fprintf(t->fd_spike_out, "%f, %d, %f\n", time + localtime, dst, t->voltage);
	}else{
	}

	localtime -= log(1.0 - drand48()) / (freq * t->neurons); // possion distributed spike timing
	dst = 1 + rand() % (t->neurons - 1); // random neuron (except reference neuron 0)
      }

      // increase time (twice because of the silence period after each noise period)
      time = (i * t->fs + j) * 2.0 * t->epochDuration;

      // increase frequency
      freq *= t->frd;
      freq += t->fad;
    } 
  }

  // close fd because fscanf sucks
  fclose(t->fd_spike_out);

  return NULL;
}