//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
//
// This generator fires gammas of fixed momentum direction (along z), from
// a random position within a square in the x-y plane at set z. The dimensions
// of the square are the same as those defining the side of the scintillator
// facing the source.
//
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......


#include "myDetPrimaryGeneratorAction.hh"
#include "myDetDetectorConstruction.hh"

#include "G4Event.hh"
#include "G4ParticleGun.hh"
#include "G4ParticleTable.hh"
#include "G4ParticleDefinition.hh"
#include "globals.hh"

#include "Randomize.hh"
#include <iostream>
#include <fstream>  // C++ files
#include <iomanip>

using namespace std;
extern ofstream fp;
extern G4double halfGasTube_length;  //declared in myDet.cc
extern G4double outerRadius;  //declared in myDet.cc

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

myDetPrimaryGeneratorAction::myDetPrimaryGeneratorAction(
                                               myDetDetectorConstruction* myDC)
:myDetector(myDC)
{
  G4int n_particle = 1;
  particleGun = new G4ParticleGun(n_particle);


  G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
  
  G4String particleName;
  particleGun->SetParticleDefinition(particleTable->FindParticle(particleName="neutron"));
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

myDetPrimaryGeneratorAction::~myDetPrimaryGeneratorAction()
{
  delete particleGun;
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void myDetPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent)
{ 

  G4float x, y, z, vx, vy, vz;
  G4ThreeVector v= G4ThreeVector(0.0,0.0,1.0);
  G4ThreeVector pos = GetRandomPosition();
 
  vx=v.getX();
  vy=v.getY();
  vz=v.getZ();
  x=pos.getX();
  y=pos.getY();
  z=pos.getZ();
  

  fp << showpoint << fixed << right;
  fp << setw(12) << setprecision(7) << vx << " "
     << setw(12) << setprecision(7) << vy << " "
     << setw(12) << setprecision(7) << vz << " "
     << setw(12) << setprecision(7) << x/cm << " "
     << setw(12) << setprecision(7) << y/cm << " "
     << setw(12) << setprecision(7) << z/cm << endl;
  particleGun->SetParticleEnergy(0.025*eV);
  particleGun->SetParticlePosition(pos);
  particleGun->SetParticleMomentumDirection(v);
  particleGun->GeneratePrimaryVertex(anEvent);
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

G4ThreeVector myDetPrimaryGeneratorAction::GetRandomPosition()
{

  G4ThreeVector retval;
  G4double rand;
  G4double xpos, ypos, zpos;


  rand = G4UniformRand();
  xpos = outerRadius * (2 * rand -1);

  rand = G4UniformRand();
  ypos = halfGasTube_length * (2 * rand -1);

  zpos = -3.*cm;

  retval.setX(xpos);
  retval.setY(ypos);
  retval.setZ(zpos);

  return retval;
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......