Correction in EMFCUT SDUM=PROD-CUT
[u/mrichter/AliRoot.git] / TFluka / TFlukaConfigOption.cxx
index 7c39dff54f668fd31818b112afd7b7cfe7ef3a7f..5f2fb095f909761ba9292ed3070792132bfdf7c0 100644 (file)
@@ -97,7 +97,22 @@ void TFlukaConfigOption::WriteFlukaInputCards()
     //
     //
     if (fMedium > -1) {
-       fprintf(fgFile,"*\n*Material specific process and cut settings for #%8d \n", fMedium);
+       TFluka* fluka = (TFluka*) gMC;
+       TObjArray *matList = fluka->GetFlukaMaterials();
+       Int_t nmaterial =  matList->GetEntriesFast();
+       TGeoMaterial* material = 0;
+       for (Int_t im = 0; im < nmaterial; im++)
+       {
+           material = dynamic_cast<TGeoMaterial*> (matList->At(im));
+           Int_t idmat = material->GetIndex();
+           if (idmat == fMedium) break;            
+       }
+       
+       
+//
+// Check if global option
+
+       fprintf(fgFile,"*\n*Material specific process and cut settings for #%8d %s\n", fMedium, material->GetName());
        fCMatMin = fMedium;
        fCMatMax = fMedium;
     } else {
@@ -382,24 +397,25 @@ void TFlukaConfigOption::ProcessCKOV()
                    cerenkovProp->GetMinimumWavelength(), cerenkovProp->GetMaximumWavelength(), cerenkovProp->GetMaximumWavelength(), 
                    Float_t(idmat), Float_t(idmat), 0.0);
            
-           if (cerenkovProp->IsMetal()) {
-               fprintf(fgFile, "OPT-PROP  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fMETAL\n", -100., -100., -100., 
-                       Float_t(idmat), Float_t(idmat), 0.0);
-           } else {
-               fprintf(fgFile, "OPT-PROP  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f\n", -100., -100., -100., 
-                       Float_t(idmat), Float_t(idmat), 0.0);
-           }
-           
+
+           fprintf(fgFile, "OPT-PROP  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f\n", -100., -100., -100., 
+                   Float_t(idmat), Float_t(idmat), 0.0);
            
            for (Int_t j = 0; j < 3; j++) {
                fprintf(fgFile, "OPT-PROP  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f&\n", -100., -100., -100., 
                        Float_t(idmat), Float_t(idmat), 0.0);
            }
-           // Photon detection efficiency user defined
-           
+
+
+           // Photon detection efficiency user defined     
            if (cerenkovProp->IsSensitive())
                fprintf(fgFile, "OPT-PROP  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fSENSITIV\n", -100., -100., -100., 
                        Float_t(idmat), Float_t(idmat), 0.0);
+           // Material has a reflective surface
+           if (cerenkovProp->IsMetal())
+               fprintf(fgFile, "OPT-PROP  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fMETAL\n", -100., -100., -100., 
+                       Float_t(idmat), Float_t(idmat), 0.0);
+
        } else {
            fprintf(fgFile,"OPT-PROD  %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fCERE-OFF\n",0., 0., 0., fCMatMin, fCMatMax, 1.);
        }
@@ -479,7 +495,7 @@ void TFlukaConfigOption::ProcessLOSS()
 // Restricted energy loss fluctuations 
 //
        fprintf(fgFile,"IONFLUCT  %10.1f%10.1f%10.1f%10.1f%10.1f\n", 1., 1., stra, fCMatMin, fCMatMax);
-       fprintf(fgFile,"DELTARAY  %10.4g%10.1f%10.1f%10.1f%10.1f%10.1f\n", cutM, 0., 0., fCMatMin, fCMatMax, 1.);       
+       fprintf(fgFile,"DELTARAY  %10.4g%10.1f%10.1f%10.1f%10.1f%10.1f\n", cutM, 0., 0., fCMatMin, fCMatMax, 1.);
     } else if (fProcessFlag[kLOSS] == 4) {
 //
 // No fluctuations
@@ -504,6 +520,7 @@ void TFlukaConfigOption::ProcessCUTGAM()
     if (fMedium == -1) {
        fprintf(fgFile,"EMFCUT    %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n", 
                0., fCutValue[kCUTGAM], 0., 0., Float_t(fgGeom->NofVolumes()), 1.);
+
     } else {
        Int_t nreg, *reglist;
        Float_t ireg;
@@ -514,6 +531,8 @@ void TFlukaConfigOption::ProcessCUTGAM()
            fprintf(fgFile,"EMFCUT    %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n", 0.,fCutValue[kCUTGAM], 0., ireg, ireg, 1.);
        }
     }
+    fprintf(fgFile,"EMFCUT    %10.4g%10.4g%10.1f%10.1f%10.1f%10.1fPROD-CUT\n", 
+           0., fCutValue[kCUTGAM], 0., fCMatMin, fCMatMax, 1.);
 }
 
 void TFlukaConfigOption::ProcessCUTELE()
@@ -534,17 +553,31 @@ void TFlukaConfigOption::ProcessCUTELE()
            fprintf(fgFile,"EMFCUT    %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n", -fCutValue[kCUTELE], 0., 0., ireg, ireg, 1.);
        }
     }
+    fprintf(fgFile,"EMFCUT    %10.4g%10.4g%10.1f%10.1f%10.1f%10.1fPROD-CUT\n", 
+           -fCutValue[kCUTELE], 0., 0., fCMatMin, fCMatMax, 1.);
 }
 
 void TFlukaConfigOption::ProcessCUTNEU()
 {
     // Cut on neutral hadrons
-    fprintf(fgFile,"*\n*Cut for neutal hadrons. CUTNEU = %13.4g\n", fCutValue[kCUTNEU]);
+    fprintf(fgFile,"*\n*Cut for neutral hadrons. CUTNEU = %13.4g\n", fCutValue[kCUTNEU]);
     if (fMedium == -1) {
        Float_t cut = fCutValue[kCUTNEU];
+       //
        // 8.0 = Neutron
        // 9.0 = Antineutron
-       fprintf(fgFile,"PART-THR  %10.4g%10.1f%10.1f\n", -cut,  8.0,  9.0);
+       //
+       // If the cut is > 19.6 MeV it is assumed the low energy neutron transport is requested.
+       // In this case the cut has to coincide with the upper  limit of the first energy group.
+       //
+       Float_t neutronCut = cut;
+       if (neutronCut < 0.0196) {
+           neutronCut = 0.0196;
+           printf("Cut on neutron lower than upper limit of first energy group.\n");
+           printf("Cut reset to 19.6 MeV !\n");
+       }
+       fprintf(fgFile,"PART-THR  %10.4g%10.1f%10.1f\n", -neutronCut,  8.0,  9.0);
+       //
        // 12.0 = Kaon zero long
        fprintf(fgFile,"PART-THR  %10.4g%10.1f%10.1f\n", -cut, 12.0, 12.0);
        // 17.0 = Lambda, 18.0 = Antilambda