// Add RANDOMIZ card
fprintf(pFlukaVmcInp,"RANDOMIZ %10.1f%10.0f\n", 1., Float_t(gRandom->GetSeed()));
// User defined ion
-// if (fUserIon) fUserIon->WriteUserInputCard(pFlukaVmcInp);
+ if (fUserIons) TFlukaIon::WriteUserInputCard(pFlukaVmcInp);
// Add START and STOP card
fprintf(pFlukaVmcInp,"START %10.1f\n",fEventsPerRun);
fprintf(pFlukaVmcInp,"STOP \n");
return (a);
}
+Int_t TFlukaIon::GetIsomerNumber(Int_t pdg)
+{
+// Acording to
+// http://cepa.fnal.gov/psm/stdhep/pdg/montecarlorpp-2006.pdf
+
+ Int_t is = pdg - 1000000000;
+ is %= 10000;
+ is %= 10;
+ return (is);
+}
+
void TFlukaIon::AddIon(Int_t a, Int_t z)
{
0, 3 * z, "Ion", pdg);
}
-void TFlukaIon::AddIon(const char* name, Int_t z, Int_t a, Int_t /*q*/,
+void TFlukaIon::AddIon(const char* name, Int_t z, Int_t a, Int_t q,
Double_t /*exE*/, Double_t mass)
{
// User defined ion
TDatabasePDG *pdgDB = TDatabasePDG::Instance();
const Double_t kAu2Gev = 0.9314943228;
- Int_t pdg = GetIonPdg(z, a);
+ Int_t pdg = GetIonPdg(z, a, q);
if (pdgDB->GetParticle(pdg)) return;
if (mass = 0.) mass = Float_t(a) * kAu2Gev + 8.071e-3;
pdgDB->AddParticle(name, "User Ion", mass, kTRUE, 0, 3 * z, "Ion", pdg);
}
-void TFlukaIon::WriteUserInputCard(FILE* pFlukaVmcInp) const
+void TFlukaIon::WriteUserInputCard(FILE* pFlukaVmcInp)
{
// Write the user input card
-
- fprintf(pFlukaVmcInp,"HI-PROPE %10.1f%10.1f%10.3f\n", (Float_t)GetZ(), (Float_t)GetA(), GetExcitationEnergy());
+ // EVENTYPE 0. 0. 2. 0. 0. 0.DPMJET
+ fprintf(pFlukaVmcInp,"EVENTYPE 0. 0. 2. 0. 0. 0.DPMJET\n");
}
Double_t GetMass() const {return fMass;}
Int_t GetPdgCode() const {return GetIonPdg(fZ, fA);}
//
- void WriteUserInputCard(FILE* file) const;
+ static void WriteUserInputCard(FILE* file);
//
static void AddIon(Int_t a, Int_t z);
static void AddIon(const char* name, Int_t z, Int_t a, Int_t q,
Double_t exE, Double_t mass);
static Int_t GetIonPdg(Int_t z, Int_t a, Int_t i = 0);
- static Int_t GetZ(Int_t pdg);
- static Int_t GetA(Int_t pdg);
+ static Int_t GetZ(Int_t pdg);
+ static Int_t GetA(Int_t pdg);
+ static Int_t GetIsomerNumber(Int_t pdg);
protected:
Int_t fZ; // Z
# define soevsv soevsv_
# define dcdion dcdion_
# define setion setion_
+# define stisbm stisbm_
#else
# define source SOURCE
# define geocrs GEOCRS
# define soevsv SOEVSV
# define dcdion DCDION
# define setion SETION
+# define stisbm STISBM
#endif
-
extern "C" {
//
// Prototypes for FLUKA functions
void type_of_call soevsv();
void type_of_call dcdion(Int_t &);
void type_of_call setion(Int_t &);
+ void type_of_call stisbm(Int_t &, Int_t &, Int_t &);
/*
*----------------------------------------------------------------------*
if (ifl == -2) {
Int_t ia = TFlukaIon::GetA(pdg);
Int_t iz = TFlukaIon::GetZ(pdg);
- IOIOCM.iproa = ia;
- IOIOCM.iproz = iz;
- BEAMCM.ijhion = iz * 1000 + ia;
- BEAMCM.ijhion = 100 * BEAMCM.ijhion + 30;
- ionid = BEAMCM.ijhion;
- dcdion(ionid);
- setion(ionid);
- FLKSTK.iloflk[FLKSTK.npflka] = ionid;
+ Int_t is = TFlukaIon::GetIsomerNumber(pdg);
+ printf("Isomer %5d \n", is);
+
+ if (is == 0) {
+ IOIOCM.iproa = ia;
+ IOIOCM.iproz = iz;
+ BEAMCM.ijhion = iz * 1000 + ia;
+ BEAMCM.ijhion = 100 * BEAMCM.ijhion + 30;
+ ionid = BEAMCM.ijhion;
+ dcdion(ionid);
+ setion(ionid);
+ FLKSTK.iloflk[FLKSTK.npflka] = ionid;
+ FLKSTK.lraddc[FLKSTK.npflka] = 0;
+ } else {
+ BEAMCM.lrdbea = 1;
+ stisbm(ia, iz, is);
+ BEAMCM.ijhion = iz * 1000 + ia;
+ BEAMCM.ijhion = 100 * BEAMCM.ijhion + 30;
+ ionid = BEAMCM.ijhion;
+ dcdion(ionid);
+ setion(ionid);
+ FLKSTK.iloflk[FLKSTK.npflka] = ionid;
+ FLKSTK.lraddc[FLKSTK.npflka] = 1;
+ }
} else {
FLKSTK.iloflk[FLKSTK.npflka] = ifl;
+ FLKSTK.lraddc[FLKSTK.npflka] = 0;
ionid = ifl;
}