modified for DeltaChord and StepMax

FullSimLight/src/FSLDetectorConstruction.cc

@@ -147,8 +147,14 @@
 #include "G4BFieldIntegrationDriver.hh"
 
 #include "G4PropagatorInField.hh"
-
+#include "G4UserLimits.hh"
+#include "G4StepLimiter.hh"
+//#include "StepMax.hh"
 #include <cassert>
+#include "HistoManager.hh"
+
+#include "G4RegionStore.hh"
+#include "G4Region.hh"
 
 //
 
@@ -173,6 +179,7 @@ FSLDetectorConstruction::FSLDetectorConstruction() : fWorld(nullptr), fDetectorM
   fMinStep             = 1.0e-2;
   fField.Put(0);
   fOutputGDMLFileName = "geometry.gdml";
+fStepMax = 2;
 //fDeltaChord;
 }
 
@@ -180,6 +187,20 @@ FSLDetectorConstruction::~FSLDetectorConstruction()
 {
   delete fDetectorMessenger;
 }
+/*/////////////////////
+FSLDetectorConstruction::FSLDetectorConstruction() : fWorld(nullptr), fDetectorMessenger(nullptr)
+{
+    // ... other member variable initializations
+
+    // Create the HistoManager and configure it
+    HistoManager* histoManager = HistoManager::GetPointer();
+    histoManager->SetFileName("y.root"); // Set the output file name
+    histoManager->SetVerbose(1); // Set verbosity level as needed
+
+}
+
+*/////////////////////
+
 
 G4VPhysicalVolume *FSLDetectorConstruction::Construct()
 {
@@ -195,7 +216,6 @@ G4VPhysicalVolume *FSLDetectorConstruction::Construct()
         if (!factory) {
             std::cout<<"Error!Cannot load geometry from factory. Exiting!"<<std::endl;
             exit(0);
-
         }
 
         GeoIntrusivePtr<GeoPhysVol> worldTmp{createGeoWorld()};
@@ -286,6 +306,68 @@ G4VPhysicalVolume *FSLDetectorConstruction::Construct()
         fWorld = physWorld;
         fWorld->GetLogicalVolume()->SetVisAttributes(G4VisAttributes::GetInvisible());
 
+//G4LogicalVolumeStore* logVolStore = G4LogicalVolumeStore::GetInstance();
+   
+/*  G4RegionStore* regionStore = G4RegionStore::GetInstance();
+  G4out <<"list of region if it worked:" << G4endl;
+    G4RegionStore* rv = (*regionStore)[i];
+    G4String regionName = rv->GetName();
+
+    // Print the logical volume's name and information
+    G4cout << "Logical Volume " << i << ": " << regionName << G4endl;    */
+/*
+G4cout << "List of all logical volumes:" << G4endl;
+for (size_t i = 0; i < logVolStore->size(); ++i) {
+    G4LogicalVolume* lv = (*logVolStore)[i];
+    G4String volumeName = lv->GetName();
+
+    // Print the logical volume's name and information
+    G4cout << "Logical Volume " << i << ": " << volumeName << G4endl;
+    // You can print additional information about the logical volume here
+
+    // If you need to access other properties of the logical volume, you can do so here
+    // For example, you can get the material or other properties:
+  //  G4Material* material = lv->GetMaterial();
+  //  G4cout << "   Material: " << material->GetName() << G4endl;
+
+    // If you need to access the daughter volumes, you can iterate through them as well:
+    for (size_t j = 0; j < lv->GetNoDaughters(); ++j) {
+        G4VPhysicalVolume* daughterVolume = lv->GetDaughter(j);
+        G4String daughterName = daughterVolume->GetName();
+        G4cout << "   Daughter Volume " << j << ": " << daughterName << G4endl;
+        // Print additional information about the daughter volume if needed
+    }
+}
+*/
+// ... (rest of your code)
+// Include necessary Geant4 headers
+// ...
+
+/*    
+    if (regionStore) {
+        G4cout << "List of defined regions:" << G4endl;
+        G4RegionStore::iterator it;
+        for (it = regionStore->begin(); it != regionStore->end(); ++it) {
+            G4Region* region = *it;
+            G4cout << "Region Name: " << region->GetName() << G4endl;
+            
+            // Access region properties 
+            G4bool isEnabled = region->IsEnabled();
+            G4bool isProductionCutSet = region->IsProductionCutSet();
+            G4double productionCut = region->GetProductionCut();
+            
+            G4cout << "   Is Enabled: " << isEnabled << G4endl;
+            G4cout << "   Is Production Cut Set: " << isProductionCutSet << G4endl;
+            G4cout << "   Production Cut: " << productionCut << G4endl;
+        }
+    } else {
+        G4cout << "No regions defined." << G4endl;
+    }
+}
+
+// ...
+*/        ///////////////////////
+
         if (fWorld == 0) {
             G4ExceptionDescription ed;
             ed << "World volume not set properly check your setup selection criteria or GDML input!" << G4endl;
@@ -302,8 +384,21 @@ G4VPhysicalVolume *FSLDetectorConstruction::Construct()
 //            std::cout<<"G4LogicalVolume store:: "<<lv->GetName()<<std::endl;
 //
 //        }
-    }
+ 
+/*     ///////////////STEPMAX/////////////////////////        
+            G4double theStepMax;  // Adjust maxStep as needed
+      G4LogicalVolumeStore *logVolStor=G4LogicalVolumeStore::GetInstance();
+        for (auto lv : *logVolStor) {				
+            G4UserLimits* userLimits = new G4UserLimits(theStepMax);
+            lv->SetUserLimits(userLimits);
+        }
 
+        // Print the value of maxStepLength (fStepMax)
+       G4cout << "StepMax (fStepMax) set to: " << theStepMax << std::endl;
+  
+
+ */
+}
     else if (fGeometryFileName.contains(".gdml")){
         G4cout << "Building the detector from the GDML file: "<<fGeometryFileName<<G4endl;
         //fParser.SetOverlapCheck(true);
@@ -426,13 +521,23 @@ void FSLDetectorConstruction::ConstructSDandField()
       G4LogicalVolumeStore *logVolStore=G4LogicalVolumeStore::GetInstance();
       for (auto lv: *logVolStore) {
 	for (auto v=plugin->beginLogVolumeNames();v!=plugin->endLogVolumeNames();v++) {
-	  if (*v==lv->GetName()) {
+	 // if (*v==lv->GetName()) {
 	    std::cout << "Sensitive Detector " << sensitiveDetector->GetName() << " attached to log vol " << lv->GetName() << " active = " << sensitiveDetector->isActive() << std::endl;
-	    lv->SetSensitiveDetector(sensitiveDetector);
-	  }
-	}
+	/*    lv->SetSensitiveDetector(sensitiveDetector);
+	    /////////////
+	    G4double theStepMax = fStepMax;  // Adjust maxStep as needed
+
+     //   for (auto lv : *logVolStore) {
+            G4UserLimits* userLimits = new G4UserLimits(theStepMax);
+            lv->SetUserLimits(userLimits);
+        }
+
+        // Print the value of maxStepLength (fStepMax)
+       std::G4cout << "StepMax (fStepMax) set to: " << fStepMax << std::endl;		
+	  }*/
+	} 
       }
-    }
+   }
 
 
 if (fFieldConstant && std::abs(fFieldValue) > 0.0) {
@@ -443,8 +548,12 @@ if (fFieldConstant && std::abs(fFieldValue) > 0.0) {
     G4FieldManager *fieldMgr = G4TransportationManager::GetTransportationManager()->GetFieldManager();
     fieldMgr->SetDetectorField(uniformMagField);
     fieldMgr->CreateChordFinder(uniformMagField);
-
-//fieldMgr->GetChordFinder()->SetDeltaChord( fDeltaChord);      
+////////STEPMAX/////
+//fieldMgr->SetMaximumTrackLength(fStepMax);  // Set the maximum step length for particle tracking
+// Set the maxStep parameter for particle tracking
+   // fieldMgr->SetMaximumTrackLength(5 * cm);  // Adjust maxStep as needed
+    // 
+    fieldMgr->GetChordFinder()->SetDeltaChord( fDeltaChord);      
 ///////////////////////////////////////////DeltaChord(mm)  Parameter Set////////////////////
 //fieldMgr->GetChordFinder()->SetDeltaChord( 0.001);   //DeltaChord(mm)  Parameter Set 1
 //fieldMgr->GetChordFinder()->SetDeltaChord( 0.1);   //DeltaChord(mm)  Parameter Set 2
@@ -471,6 +580,28 @@ else {
     }
 }
 
+//
+// Define the maximum step length
+//G4double theStepMax = 1.0 * mm; // Set your desired maximum step length (e.g., 1 mm)
+G4double theStepMax;
+// Get the logical volume store
+G4LogicalVolumeStore* logVolStore = G4LogicalVolumeStore::GetInstance();
+
+// Loop through logical volumes and set the step limits
+for (auto lv : *logVolStore) {
+    G4UserLimits* userLimits = lv->GetUserLimits();
+    if (!userLimits) {
+        userLimits = new G4UserLimits(theStepMax);
+        lv->SetUserLimits(userLimits);
+    } else {
+        // If user limits already exist, update the maximum step length
+        userLimits->SetMaxAllowedStep(theStepMax);
+    }
+}
+
+// Print the value of maxStepLength (theStepMax)
+G4cout << "StepMax (theStepMax) set to: " << theStepMax << std::endl;
+   
 
 // This is thread-local
 //To print the values of parameters