How to Simulate charged-particle tracks using Heed. i have tryied many time but not getting result,

import ROOT
import gmsh
import ctypes
import sys
import math
import numpy as np
import matplotlib.pyplot as plt
gmsh.initialize()

Set relevant LEM parameters.

LEM thickness in cm

lem_th = 0.04;

copper thickness

lem_cpth = 0.0035;
#LEM pitch in cm
lem_pitch = 0.07;

X-width of drift simulation will cover between +/- axis_x

axis_x = 0.1;

Y-width of drift simulation will cover between +/- axis_y

axis_y = 0.1;
axis_z = 0.25 + lem_th / 2 + lem_cpth;

Define the medium.

ROOT.gSystem.Load(“/home/dileshwar/garfieldpp/install/lib/libmagboltz.so.11”)
ROOT.gSystem.Load(“/home/dileshwar/garfieldpp/install/lib/libGarfield.so”)
gas = ROOT.Garfield.MediumMagboltz()
gas.SetTemperature(293.15) #Set the temperature [K]
gas.SetPressure(740.)
gas.EnableDrift()
gas.SetComposition(“Ar”, 70., “co2”, 30.)#Specify the gas mixture (Ar/CO2 70:30)s.
#gas.LoadIonMobility(‘/share/Garfield/Data/IonMobility_Ar+_Ar.txt’)
gas.Initialise(True)

Define the gas medium using Magboltz.

gas = ROOT.Garfield.MediumMagboltz()
gas.SetTemperature(293.15) # in Kelvin
gas.SetPressure(740.) # in Torr
gas.EnableDrift() # Enable electron drift in the medium
gas.SetComposition(“Ar”, 70., “co2”, 30.) # Set gas mixture composition (70% Ar, 30% CO2)

Read in the 3D field map.

elm = ROOT.Garfield.ComponentElmer(
“gemcell/mesh.header”, “gemcell/mesh.elements”, “gemcell/mesh.nodes”,
“gemcell/dielectrics.dat”, “gemcell/gemcell.result”, “cm”
)
elm.EnablePeriodicityX() # Enable periodicity in X
elm.EnableMirrorPeriodicityY() # Enable mirror periodicity in Y
elm.SetMedium(0, gas) # Assign the medium to region 0

Create a Sensor object.

sensor = ROOT.Garfield.Sensor()
sensor.AddComponent(elm)
sensor.SetArea(-axis_x, -axis_y, -axis_z, axis_x, axis_y, axis_z)
track = ROOT.Garfield.TrackHeed()
track.SetSensor(sensor) # Set the sensor for Heed
track.SetParticle(“muon”) # You can change this to the desired particle type (e.g., “electron”, “muon”)
track.SetEnergy(180.e9) # Set the particle energy (in eV). Adjust as needed.

drift = ROOT.Garfield.DriftLineRKF()
drift.SetGainFluctuationsPolya(0., 20000.)
drift.SetSensor(sensor)
driftView = ROOT.Garfield.ViewDrift()
cD = ROOT.TCanvas(“cD”, “”, 600, 600)
driftView.SetCanvas(cD)

drift.EnablePlotting(driftView)
track.EnablePlotting(driftView)
zi = 0.5 * lem_th + lem_cpth + 0.1 # Starting z position
ri = (lem_pitch / 2) * np.random.uniform() # Random radius within the LEM pitch
thetai = np.random.uniform() * 2 * np.pi # Random angle
xi = ri * np.cos(thetai) # x position based on random angle
yi = ri * np.sin(thetai)
track.NewTrack(xi, yi, zi, 0., 0., 0., 0.) # Assume particle moving in -z direction

Loop over the clusters along the track.

for cluster in track.GetClusters():
# Loop over the electrons in the cluster.
for electron in cluster.electrons:
drift.DriftElectron(electron.x , electron.y, electron.z, electron.t)

Draw the canvas to visualize drift lines and tracs

#cD.Draw()

Visualize the drift lines and field using ViewFEMesh.

cGeom = ROOT.TCanvas(“cD”, “”)
viewMesh = ROOT.Garfield.ViewFEMesh()
viewMesh.SetArea(-axis_x, -axis_z, -axis_y, axis_x, axis_z, axis_y)
viewMesh.SetCanvas(cGeom)
viewMesh.SetComponent(elm)
viewMesh.SetPlane(0, -1, 0, 0, 0, 0) # Set plane for plotting
viewMesh.SetFillMesh(True)
viewMesh.SetColor(1, ROOT.kGray)
viewMesh.SetColor(2, ROOT.kYellow + 3)
viewMesh.SetColor(3, ROOT.kYellow + 3)
viewMesh.EnableAxes()
viewMesh.SetXaxisTitle(“x (cm)”)
viewMesh.SetYaxisTitle(“y (cm)”)
viewMesh.SetViewDrift(driftView)
viewMesh.Plot()
cGeom.Draw()

Hi,
did you try one of the examples on the repository that use Heed? For instance this one?

Or this one:

Or this one:

1 Like

Thank you for your quick reply! i haven’t use examples on the repository that use Heed. Can you tell me what is wrong in code, please give me suggestions.

Can you tell us what is the error or the non-desired/expected result you obtained?
Without your input files (“gemcell/mesh.header”, “gemcell/mesh.elements”, “gemcell/mesh.nodes”, “gemcell/dielectrics.dat”, “gemcell/gemcell.result”) we can not run your example.
Kind regards
Piet

Indeed, you didn’t really gave us any hint as to what error you got. Having said that, looking at your code I believe the issue is that you didn’t load a gas file. DriftLineRKF (and AvalancheMC) need pre-computed tables of transport coefficients like drift velocity, diffusion coefficients, etc. to simulate electron drift lines.

i am encountering with error
MediumMagboltz::SetComposition: Ar/CO2 (70/30)
ComponentElmer::Initialise:
Read 97248 nodes and 67732 elements from file gemcell/mesh.header.
Read potentials from file gemcell/gemcell.result.
Set material 0 of 4 to eps 1.
Set material 1 of 4 to eps 3.23.
Set material 2 of 4 to eps 1e+10.
Set material 3 of 4 to eps 1e+10.
ComponentElmer::Initialise: Finished.
ComponentElmer::Prepare:
Caching the bounding boxes of all elements… done.
Initialized tetrahedral tree.
DriftLineRKF::SetGainFluctuationsPolya: Mean avalanche size set to 20000.

DriftLineRKF::GetVelocity:
Cannot retrieve drift velocity at (-0.008794, 0.010352, 0.135514).
DriftLineRKF::DriftLine:
Cannot retrieve drift velocity at initial position (-0.008794, 0.010352, 0.135514).
DriftLineRKF::GetVelocity:
Cannot retrieve drift velocity at (-0.008634, 0.008446, 0.137157).
DriftLineRKF::DriftLine:
Cannot retrieve drift velocity at initial position (-0.008634, 0.008446, 0.137157).
DriftLineRKF::GetVelocity:
Cannot retrieve drift velocity at (-0.006294, -0.020242, 0.161458).
DriftLineRKF::DriftLine:
Cannot retrieve drift velocity at initial position (-0.006294, -0.020242, 0.161458).
DriftLineRKF::GetVelocity:
Cannot retrieve drift velocity at (-0.004516, -0.042723, 0.180503).
DriftLineRKF::DriftLine:
Cannot retrieve drift velocity at initial position (-0.004516, -0.042723, 0.180503).
DriftLineRKF::GetVelocity:
Cannot retrieve drift velocity at (-0.003966, -0.049413, 0.186216).
DriftLineRKF::DriftLine:
Cannot retrieve drift velocity at initial position (-0.003966, -0.049413, 0.186216).

Indeed, as @hschindl pointed out, this is a typical error message you get when you do not have set the correct gas medium for DriftLineRKF. He also pointed you to the examples, because for instance here you can see that before using DriftLineRKF one loads a gasfile.

Thanks for support & understanding :pray: