//
// Program: neutron.C
//
// Description: Simulazione di un flusso di neutroni attraverso la materia
//
// Author: A. Garfagnini
//
// Date: 30 November 2013
//
#include "TApplication.h"     // SA version
#include "TROOT.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TH2F.h"
#include "TGraph.h"
#include "TLine.h"
#include "TLegend.h"
#include "TBox.h"
#include "TArrow.h"

#include "boost/random.hpp"


#include <iostream>
#include <cmath>
#include <fstream>
#include <unistd.h>
#include <libgen.h>
#include <stdlib.h>

using namespace std;

void usage(const char * pname)
{
    cerr << "Usage: " << pname << " [-g event]" << endl;
}

void neutron()
{
  gROOT->SetStyle("Plain");
  gStyle->SetOptStat(0);
  gStyle->SetOptTitle(kFALSE);
  gStyle->SetLabelOffset(0.01,"x");
  gStyle->SetLabelOffset(0.005,"y");
  gStyle->SetTitleOffset(1.2,"y");
  gStyle->SetTitleOffset(1.2,"x");
}


int main(int argc, char * argv[])
{
   int opt;
   int display_event = 0;

    // Handle command line options
    while((opt = getopt(argc, argv, "g:h?")) != -1) {

        switch (opt) {

        case 'g':
            display_event = strtol(optarg, NULL, 0);
            break;

        case 'h':
        case '?':
        default:
            usage( basename( argv[0] ) );
            return 1;
        }
    }

    // Handle user data input
    cout << "Random generator seed: ";
    int seed;
    cin >> seed;

    //boost::uniform_real<> uni(-M_PI, M_PI);
    boost::uniform_real<> uni(-1.0, 1.0);
    boost::mt19937 rng;
    rng.seed(seed);
    boost::variate_generator < boost::mt19937 &,
                               boost::uniform_real<> > uni_one(rng, uni);


    // Startup ROOT
    TApplication myapp("myapp", &argc, argv);

    gROOT->Reset();
    neutron();

    TCanvas * c3 = new TCanvas("c3","test");
    gPad->SetTopMargin(0.01);
    gPad->SetBottomMargin(0.12);
    gPad->SetRightMargin(0.01);
    gPad->SetLeftMargin(0.12);

    TH2F * hs = new TH2F("hs", "", 50, -1., 6., 50, -8., 8.0);
    hs->Draw();
    hs->GetXaxis()->SetTitle("ntry");
    hs->GetXaxis()->SetLabelSize(0.04);
    hs->GetYaxis()->SetTitle("result");

    TBox * hb = new TBox(0., -8.0, 5.0, 8.0);
    hb->Draw();

    double x[19] = {0.0};
    double y[19] = {0.0};

    // Segui 18 interazioni successive
    x[0] = 0.0;
    y[0] = 0.0;
    x[1] = 1.0;
    y[1] = 0.0;
    TArrow * ha = new TArrow(x[0], y[0], x[1], y[1]);
    ha->SetLineColor(4);
    ha->Draw();
    for (int j=2; j<18; j++) {

        double u = uni_one();
        double theta = acos(u);
        x[j] = x[j-1] + cos(theta);
        y[j] = y[j-1] + sin(theta);

        TArrow * ha = new TArrow(x[j-1], y[j-1], x[j], y[j]);
        ha->SetLineColor(4);
        ha->Draw();

        if (x[j] <= 0.0) {
            cout << "neutrone back-scattered" << endl;
            break;
        }

        if (x[j] >= 5.0) {
            cout << "neutrone uscito" << endl;
            break;
        }
    }

    c3->Update();

    myapp.Run();

    return 0;
}