// @(#)CTDR/Electron:$Name:  $:$Id: CTDR_Electron,cpp,v 2.0 2006/11/24 16:10:45 brun Exp $
// Author: Mathieu Benoit   24/11/06

/*************************************************************************
 * . Class CTDR_Interaction				                 *
 *   Written by : Mathieu Benoit                                         *
 *   m.benoit@umontreal.ca						 *
 *************************************************************************/
#include <cstdlib>
#include <iostream>
#include <fstream>
#include <math.h>
#include <stdlib.h>
#include <time.h>
#include <stdio.h>
#include <string>
#include "CTDR_Interaction.h"
using namespace std;

//////////////////////////////////////////////////////////////////////////
// Source file : class CTDR_Interaction					//	      
// This class represent a cloud of charge element of a gamma 		//
// CTDR_Interaction in CZT.It contains all the function to create 	//
// and simulate the behavior of an carrier cloud generated by a 	//
// gamma event. The cloud consiste of CTDR_Electron objects,		//
// wich all evolves in time from their generation to their collection. 	//		      
// 									//
// Constant to set in source before compilation				//
// "Nelec" represent the number of charge element (class CTDR_Electrons)// 
// used to reproduce the real charge cloud generated by a gamma 	//
// CTDR_Interaction. The higher the number, he more accurately the 	//
// behavior is simulated. However , higher number of charge element also//
// mean longer computation time. CTDR_Fields and weighting potentials 	//
// must be extracted from solution files obtained with poisson equation //
// solvers for each charge at each time step.				//
// To simulate signals and to focus on one event or one multi-hit, 	//
// a high number of charge is reccomended. But for simulation of a large// 
// number of events, for example to obtain scatter plots, a small number// 
// of charge works fine and give much more reasonable. As computer 	//
// become faster, these type of simulation will gain much precision as 	//
// they require a lot of brutal force.These simulation can be performed //
// on a single-processor , but the time needed to obtain readable data 	//
// is quite long,(days or weeks). These simulations can easily be 	//
// performed in parallel on many processor.The average parallel		//
// computing facility these days could obtain good results quite fast 	//
// for simulation with a high number of charge elements. This code 	//
// can easily scale to the demand and to the ressource at your 		//
// disposition to give  good results in a maximum of situations.	//
//									//
// zmax is the maximum height of the CZT block 				//
//////////////////////////////////////////////////////////////////////////
const int Nelec=100;
const double  zmax=0.0075;
const double pi=M_PI;

ClassImp(CTDR_Interaction)
//______________________________________________________________________________
CTDR_Interaction::CTDR_Interaction(double x, double y, double z,double Energie,double sigma,double t0)
{
//Constructor for the CTDR_Interaction , each CTDR_Interaction built represent 
//one CTDR_Interaction site with deposited energy "Energie", and the cloud 
//begins its simulation with a radius of sigma. It has been found that 
//Electron-hole separation process is fast and can be neglected.. , 

	double xtemp; // Centroid position in X
	double ytemp; // Centroid position in Y
	double ztemp;// Centroid position in Z
	double Rtemp,R2temp; // Random variable
	energie=Energie;  // Interaction's energy
	N=Nelec; //The number of charge elements used for the simulation
	t=t0;

	// The charge elements are distibuted according to a gaussian of sigma in each dimension
	// Around the CTDR_Interaction site position
	for(int i=0;i<Nelec;i++){
		Rtemp=(float)rand()/RAND_MAX;
        	R2temp=(float)rand()/RAND_MAX;
		xtemp=x+sigma*sin(2*pi*Rtemp)*sqrt(-2*log(R2temp));
		Rtemp=(float)rand()/RAND_MAX;
        	R2temp=(float)rand()/RAND_MAX;
		ytemp=y+sigma*sin(2*pi*Rtemp)*sqrt(-2*log(R2temp));
		Rtemp=(float)rand()/RAND_MAX;
        	R2temp=(float)rand()/RAND_MAX;
		ztemp=z+sigma*sin(2*pi*Rtemp)*sqrt(-2*log(R2temp));
		//We set the mobility, 0.1 for CTDR_Electron in CZT, and the amount of real charge 
		// attributed to each charge elements
		charges[i].SetValeurInitial
			(xtemp,ytemp,ztemp,0.1,energie/(Nelec*4.4));
		};
};

//______________________________________________________________________________
bool CTDR_Interaction::iteration(double dt,CTDR_Field &A)
{
//This function perform a time step for all the charge element in the cloud 

       bool I=0;
		for(int i=0;i<Nelec;i++){
		if(charges[i].GetZ()<zmax){ //If charge have not reached the anode

		//Charges cannot go lower than the cathode
		if(charges[i].GetZ()<0){charges[i].SetPosition(charges[i].GetX(),charges[i].GetY(),0.0000001);};
		//We do a time step for each charge element
		charges[i].Iterateur(dt,int(energie/4.4),rmsX,rmsY,rmsZ,A);
			}
		else{
			//We set the charge as collected
			charges[i].Collecte(zmax);


		        };
	//If no charge are free, I is 1, else its 0
        if(charges[i].GetEcollecte()==Libre){I=1;};
		};

		t+=dt; //time advance

		//If all charge are collected , we return the signal to stop for this CTDR_Interaction to the global simulation
		if(I==0){return 1;}
		else {return 0;};
		//cout << t << endl;

};

//______________________________________________________________________________
void CTDR_Interaction::ComputeRms()
{
//This function computes the RMS in x,y,z of the charge cloud, the centroid is also calculated
//These data are passed at each time step the the charge elements to set-up their random walk 

	double carrex=0,carrey=0,carrez=0,moyx=0,moyy=0,moyz=0;
	for(int i=0;i<Nelec;i++) {
		carrex+=pow(charges[i].GetX(),2);
		carrey+=pow(charges[i].GetY(),2);
		carrez+=pow(charges[i].GetZ(),2);
		moyx+=charges[i].GetX();
		moyy+=charges[i].GetY();
		moyz+=charges[i].GetZ();
		};
	rmsX=sqrt(fabs(carrex/Nelec - pow(moyx/Nelec,2)));
	rmsY=sqrt(fabs(carrey/Nelec - pow(moyy/Nelec,2)));
	rmsZ=sqrt(fabs(carrez/Nelec - pow(moyz/Nelec,2)));
	MoyZ=moyz;
	//cout << t << ' ' << rmsX << ' ' << rmsY << ' ' << rmsZ << endl ;
};

//______________________________________________________________________________
double CTDR_Interaction::GetRmsX()
{
//Acces method to RMS

	return rmsX;
};

//______________________________________________________________________________
double CTDR_Interaction::GetRmsY()
{
//Acces method to RMS
	return rmsY;
};

//______________________________________________________________________________
double CTDR_Interaction::GetRmsZ()
{
//Acces method to RMS
	return rmsZ;
};

//______________________________________________________________________________
double CTDR_Interaction::GetMoyZ()
{
//Acces method to centroid position in Z
	return MoyZ/Nelec;
};

//______________________________________________________________________________
double CTDR_Interaction::GetQx()
{
//Acces method to centroid position in Z

	return Q0x;
};

//______________________________________________________________________________
double CTDR_Interaction::GetQy()
{
//Acces method to centroid position in Z

	return Q0y;
};

//______________________________________________________________________________
void CTDR_Interaction::SetQx(double Q)
{
//Modifier function of the initial weight of the charge cloud. Data used in signal computation

	Q0x=Q;
};

//______________________________________________________________________________
void CTDR_Interaction::SetQy(double Q)
{
//Modifier function of the initial weight of the charge cloud. Data used in signal computation

	Q0y=Q;
};

//______________________________________________________________________________
int CTDR_Interaction::GetN()
{
// Access method for Nelec

return N;
};

//______________________________________________________________________________
CTDR_Electron CTDR_Interaction::GetCharge(int i)
{
// Access method for Electrons object
	return charges[i];
};

//______________________________________________________________________________
double CTDR_Interaction::Gett()
{
// Access method for time
	return t;
};

//______________________________________________________________________________
void CTDR_Interaction::WritePositions(ofstream &sortie)
{
// This function write the position of the charge elements to a file


	for(int i=0;i<Nelec;i++){
		sortie << charges[i].GetX() << ' ' << charges[i].GetY() << ' '
		<< charges[i].GetZ() << ' ' << charges[i].GetQ() << endl;
		};
};