G4AblaDataDefs.hh

Go to the documentation of this file.

//
// ********************************************************************
// * License and Disclaimer                                           *
// *                                                                  *
// * The  Geant4 software  is  copyright of the Copyright Holders  of *
// * the Geant4 Collaboration.  It is provided  under  the terms  and *
// * conditions of the Geant4 Software License,  included in the file *
// * LICENSE and available at  http://cern.ch/geant4/license .  These *
// * include a list of copyright holders.                             *
// *                                                                  *
// * Neither the authors of this software system, nor their employing *
// * institutes,nor the agencies providing financial support for this *
// * work  make  any representation or  warranty, express or implied, *
// * regarding  this  software system or assume any liability for its *
// * use.  Please see the license in the file  LICENSE  and URL above *
// * for the full disclaimer and the limitation of liability.         *
// *                                                                  *
// * This  code  implementation is the result of  the  scientific and *
// * technical work of the GEANT4 collaboration.                      *
// * By using,  copying,  modifying or  distributing the software (or *
// * any work based  on the software)  you  agree  to acknowledge its *
// * use  in  resulting  scientific  publications,  and indicate your *
// * acceptance of all terms of the Geant4 Software license.          *
// ********************************************************************
//
// ABLAXX statistical de-excitation model
// Pekka Kaitaniemi, HIP (translation)
// Christelle Schmidt, IPNL (fission code)
// Davide Mancusi, CEA (contact person INCL/ABLA)
// Aatos Heikkinen, HIP (project coordination)
//
#define ABLAXX_IN_GEANT4_MODE 1

#include "globals.hh"

// Data structures needed by ABLA evaporation code.

#ifndef G4AblaDataDefs_hh
#define G4AblaDataDefs_hh 1

#ifdef ABLAXX_IN_GEANT4_MODE
#include "globals.hh"
#else
#include "G4INCLGeant4Compat.hh"
#endif

#include <cmath>

// ABLA

class G4Nevent {
public:
  G4Nevent() {};
  ~G4Nevent() {};
  
  G4int ii;
};

// ABLA
#define PACESIZEROWS 500
#define PACESIZECOLS 500

class G4Pace {

public:
  G4Pace() {};

  ~G4Pace() {};
  
  G4double dm[PACESIZEROWS][PACESIZECOLS];
};

#define EC2SUBROWS 154
#define EC2SUBCOLS 99

class G4Ec2sub {
public:
  G4Ec2sub() {};

  ~G4Ec2sub() {};

  G4double ecnz[EC2SUBROWS][EC2SUBCOLS]; 

  void dump() {
    for(G4int i = 0; i < EC2SUBROWS; i++) {
      for(G4int j = 0; j < EC2SUBCOLS; j++) {
        //G4cout << ecnz[i][j] << " ";
      }
      //      G4cout << G4endl;
    }
  }
};

class G4Ald {
public:
  G4Ald()
    :av(0.0), as(0.0), ak(0.0), optafan(0.0)
  {};
  ~G4Ald() {};
  
  G4double av,as,ak,optafan;
};

#define ECLDROWS 154
#define ECLDCOLS 99

class G4Ecld {

public:
  G4Ecld() {};
  ~G4Ecld() {};

  G4double ecgnz[ECLDROWS][ECLDCOLS];

  G4double ecfnz[ECLDROWS][ECLDCOLS];

  G4double vgsld[ECLDROWS][ECLDCOLS]; 

  G4double alpha[ECLDROWS][ECLDCOLS];
};

class G4Fiss {
public:
  G4Fiss()
    :akap(0.0), bet(0.0), homega(0.0), koeff(0.0), ifis(0.0),
     optshp(0), optxfis(0), optles(0), optcol(0)
  {};
  ~G4Fiss() {};
  
  G4double akap,bet,homega,koeff,ifis;
  G4int optshp, optxfis,optles,optcol;
};

#define FBROWS 101
#define FBCOLS 161

class G4Fb {
  
public:
  G4Fb() {};
  ~G4Fb() {;}
  
  //  G4double efa[FBROWS][FBCOLS];
  G4double efa[FBCOLS][FBROWS];
};

class G4Opt {

public:
  G4Opt()
    :optemd(0), optcha(0), eefac(0.0)
  {};
  ~G4Opt() {};
  
  G4int optemd,optcha;
  G4double eefac;                                  
};

#define EENUCSIZE 2002
#define XHESIZE 50
class G4Eenuc {
public:
  G4Eenuc() {
    for(G4int i = 0; i < EENUCSIZE; ++i) {
      she[i] = 0.0;
    }
    for(G4int i = 0; i < XHESIZE; ++i) {
      for(G4int j = 0; j < EENUCSIZE; ++j) {
        xhe[i][j] = 0.0;
      }
    }
  };
  ~G4Eenuc() {};
  
  G4double she[EENUCSIZE],xhe[XHESIZE][EENUCSIZE];                                            
};

//#define VOLANTSIZE 200
#define VOLANTSIZE 301

class G4Volant {
  
public:
  G4Volant()
  {
    clear();
  }

  ~G4Volant() {};

  void clear()
  {
    for(G4int i = 0; i < VOLANTSIZE; i++) {
      copied[i] = false;
      acv[i] = 0;
      zpcv[i] = 0;
      pcv[i] = 0;
      xcv[i] = 0;
      ycv[i] = 0;
      zcv[i] = 0;
      iv = 0;
    }
  }

  G4double getTotalMass()
  {
    G4double total = 0.0;
    for(G4int i = 0; i <= iv; i++) {
      total += acv[i];
    }
    return total;
  }

  void dump()
  {
    G4double totA = 0.0, totZ = 0.0, totP = 0.0;
    //    G4cout <<"i \t ACV \t ZPCV \t PCV" << G4endl; 
    for(G4int i = 0; i <= iv; i++) {
      if(i == 0 && acv[i] != 0) {
        //      G4cout <<"G4Volant: Particle stored at index " << i << G4endl;
      }
      totA += acv[i];
      totZ += zpcv[i];
      totP += pcv[i];
      //      G4cout << "volant" << i << "\t" << acv[i] << " \t " << zpcv[i] << " \t " << pcv[i] << G4endl;
    }
    //    G4cout <<"Particle count index (iv) = " << iv << G4endl;
    //    G4cout <<"ABLA Total: A = " << totA << " Z = " << totZ <<  " momentum = " << totP << G4endl;
  }

  G4double acv[VOLANTSIZE],zpcv[VOLANTSIZE],pcv[VOLANTSIZE],xcv[VOLANTSIZE];
  G4double ycv[VOLANTSIZE],zcv[VOLANTSIZE];
  G4bool copied[VOLANTSIZE];
  G4int iv; 
};

#define VARNTPSIZE 301
class G4VarNtp {
public:
  G4VarNtp() {
    clear();
  };

  ~G4VarNtp() {};

  void clear() {
    particleIndex = 0;
    projType = 0;
    projEnergy = 0.0;
    targetA = 0;
    targetZ = 0;
    masp = 0.0; mzsp = 0.0; exsp = 0.0; mrem = 0.0;
    // To be deleted?
    spectatorA = 0;
    spectatorZ = 0;
    spectatorEx = 0.0;
    spectatorM = 0.0;
    spectatorT = 0.0;
    spectatorP1 = 0.0;
    spectatorP2 = 0.0;
    spectatorP3 = 0.0;
    massini = 0;
    mzini = 0;
    exini = 0;
    pcorem = 0;
    mcorem = 0;
    pxrem = 0;
    pyrem = 0;
    pzrem = 0;
    erecrem = 0;
    mulncasc = 0;
    mulnevap = 0;
    mulntot = 0;
    bimpact = 0.0;
    jremn = 0;
    kfis = 0;
    estfis = 0;
    izfis = 0;
    iafis = 0;
    ntrack = 0;
    needsFermiBreakup = false;
    for(G4int i = 0; i < VARNTPSIZE; i++) {
      itypcasc[i] = 0;
      avv[i] = 0;
      zvv[i] = 0;
      enerj[i] = 0.0;
      plab[i] = 0.0;
      tetlab[i] = 0.0;
      philab[i] = 0.0;
      full[i] = false;
    }
  }

  void addParticle(G4double A, G4double Z, G4double E, G4double P, G4double theta, G4double phi) {
    if(full[particleIndex]) {
      //      G4cout <<"A = " << Z << " Z = " << Z << G4endl;
    } else {
      avv[particleIndex] = (int) A;
      zvv[particleIndex] = (int) Z;
      enerj[particleIndex] = E;
      plab[particleIndex] = P;
      tetlab[particleIndex] = theta;
      philab[particleIndex] = phi;
      full[particleIndex] = true;
      ntrack = particleIndex + 1;
      particleIndex++;
    }
  }

  G4int getTotalBaryonNumber() {
    G4int baryonNumber = 0;
    for(G4int i = 0; i < ntrack; i++) {
      if(avv[i] > 0) {
        baryonNumber += avv[i];
      }
    }
    return baryonNumber;
  }

  G4double getTotalEnergy() {
    G4double energy = 0.0;
    for(G4int i = 0; i < ntrack; i++) {
      energy += std::sqrt(std::pow(plab[i], 2) + std::pow(getMass(i), 2)); // E^2 = p^2 + m^2
    }

    return energy;
  }

  G4double getTotalThreeMomentum() {
    G4double momentum = 0;
    for(G4int i = 0; i < ntrack; i++) {
      momentum += plab[i];
    }
    return momentum;
  }

  G4double getMomentumSum() {
    G4double momentum = 0;
    for(G4int i = 0; i < ntrack; i++) {
      momentum += plab[i];
    }
    return momentum;
  }

  G4double getMass(G4int particle) {
    const G4double protonMass = 938.272;
    const G4double neutronMass = 939.565;
    const G4double pionMass = 139.57;

    G4double mass = 0.0;
    if(avv[particle] ==  1 && zvv[particle] ==  1) mass = protonMass;
    if(avv[particle] ==  1 && zvv[particle] ==  0) mass = neutronMass;
    if(avv[particle] == -1)                        mass = pionMass;
    if(avv[particle] > 1)
      mass = avv[particle] * protonMass + zvv[particle] * neutronMass;
    return mass;
  }

  void dump() {
    G4int nProton = 0, nNeutron = 0;
    G4int nPiPlus = 0, nPiZero = 0, nPiMinus = 0;
    G4int nH2 = 0, nHe3 = 0, nAlpha = 0;
    G4int nFragments = 0;
    G4int nParticles = 0;
    for(G4int i = 0; i < ntrack; i++) {
      nParticles++;
      if(avv[i] ==  1 && zvv[i] ==  1) nProton++;  // Count multiplicities
      if(avv[i] ==  1 && zvv[i] ==  0) nNeutron++;
      if(avv[i] == -1 && zvv[i] ==  1) nPiPlus++;
      if(avv[i] == -1 && zvv[i] ==  0) nPiZero++;
      if(avv[i] == -1 && zvv[i] == -1) nPiMinus++;
      if(avv[i] ==  2 && zvv[i] ==  1) nH2++;
      if(avv[i] ==  3 && zvv[i] ==  2) nHe3++;
      if(avv[i] ==  4 && zvv[i] ==  2) nAlpha++;
      if(                zvv[i] >   2) nFragments++;
    }
  }

  G4int projType;

  G4double projEnergy;

  G4int targetA;

  G4int targetZ;

  G4double masp, mzsp, exsp, mrem;

  G4int spectatorA;

  G4int spectatorZ;

  G4double spectatorEx;

  G4double spectatorM;

  G4double spectatorT;

  G4double spectatorP1;

  G4double spectatorP2;

  G4double spectatorP3;

  G4double massini;

  G4double mzini;

  G4double exini;

  G4double pcorem, mcorem, pxrem, pyrem, pzrem, erecrem;

  G4int mulncasc;

  G4int mulnevap;

  G4int mulntot;

  G4double bimpact;

  G4int jremn;

  G4int kfis;

  G4double estfis;

  G4int izfis;

  G4int iafis;

  G4int ntrack;

  G4bool full[VARNTPSIZE];

  G4bool needsFermiBreakup;

  G4int itypcasc[VARNTPSIZE];

  
  G4int avv[VARNTPSIZE];

  G4int zvv[VARNTPSIZE];

  G4double enerj[VARNTPSIZE];

  G4double plab[VARNTPSIZE];

  G4double tetlab[VARNTPSIZE];

  G4double philab[VARNTPSIZE];

private:
  G4int particleIndex;
};

#endif