Isomeres are supported.
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 18 Apr 2008 12:10:27 +0000 (12:10 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 18 Apr 2008 12:10:27 +0000 (12:10 +0000)
TFluka/TFluka.cxx
TFluka/TFlukaIon.cxx
TFluka/TFlukaIon.h
TFluka/source.cxx

index eaeab8d..f384afd 100644 (file)
@@ -1330,7 +1330,7 @@ void TFluka::InitPhysics()
 // 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");
index e888a08..09af298 100644 (file)
@@ -82,6 +82,17 @@ Int_t TFlukaIon::GetA(Int_t pdg)
     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)
 {
 
@@ -95,22 +106,22 @@ 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");
 }
index 989b642..95282bd 100644 (file)
@@ -32,14 +32,15 @@ public:
     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
index 768e68c..c07ec57 100644 (file)
@@ -33,6 +33,7 @@
 # define soevsv soevsv_
 # define dcdion dcdion_
 # define setion setion_
+# define stisbm stisbm_
 #else
 # define source SOURCE
 # define geocrs GEOCRS
@@ -41,8 +42,8 @@
 # define soevsv SOEVSV
 # define dcdion DCDION
 # define setion SETION
+# define stisbm STISBM
 #endif
-
 extern "C" {
   //
   // Prototypes for FLUKA functions
@@ -54,6 +55,7 @@ extern "C" {
   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 &);
 
  /*
    *----------------------------------------------------------------------*
@@ -197,16 +199,33 @@ extern "C" {
        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;
        }