Drift lines problem with trackheed

Dear experts,

I’m currently working on protons crossing a plane parallel ionisation chamber geometry and studying the attachment. I have some problems when i try to understand the visualization of the drift lines.

  1. Why my proton (green line) stop in the chamber while the energy it loses is negligible ?
  2. Why does it look like there are 2 electrons, one collected by the electrode and the other attached, when the code give a total number of electrons of 1?
  3. In the second photo, i have a drift line of electron who stop but all the status are -7 (attachment) or -5 (lft the drift medium =? collected by the electrode) you can see the status printed in the terminal to the right.

One last question: is it possible to change the color of the positive or negative ion to distinguish them ?


Thank you for your help.

Pierre

// Make a medium
MediumMagboltz gas("o2", 21.1, "Ar", .9, "H2O", 2., "N2", 76.);
gas.SetTemperature(273.15 + 20);
gas.SetPressure(760.);
const std::string path = getenv("GARFIELD_INSTALL");
gas.LoadIonMobility(path + "/share/Garfield/Data/IonMobility_N2+_N2.txt");
gas.LoadNegativeIonMobility(path + "/share/Garfield/Data/IonMobility_O2-_air.txt");

// Thickness of the gas gap [cm] and voltage [V]
constexpr double gap = .1;
constexpr double voltage = 500.;

// Drift canvas
const bool plotDrift = true;
TCanvas* cD = nullptr;
if (plotDrift) cD = new TCanvas("cD", "", 600, 600);

// Make a component
ComponentConstant cmp;
cmp.SetArea(0., -1., -1., gap, 1., 1.);
cmp.SetMedium(&gas);
double field = voltage / double(gap);
cmp.SetElectricField(field, 0., 0.);

// Make a sensor
Sensor sensor;
sensor.AddComponent(&cmp);

// DriftView class
ViewDrift driftView;
float ViewWidth = .02;
driftView.SetArea(0., -ViewWidth, -ViewWidth, gap, ViewWidth, ViewWidth);
if (plotDrift) driftView.SetCanvas(cD);

// AvalancheMicroscopic class (for electrons)
AvalancheMicroscopic aval;
aval.SetSensor(&sensor);
aval.EnablePlotting(&driftView);
aval.EnableExcitationMarkers(true);


// AvalancheMC class (for ions)
AvalancheMC drift;
drift.SetSensor(&sensor);
drift.EnablePlotting(&driftView);

// Track class
TrackHeed track;
track.SetSensor(&sensor);
track.SetParticle("proton");
//track.SetKineticEnergy(230.e6); // 230 MeV
track.SetGamma(3); // using MIP for more visibility in the drift view
constexpr bool verbose = true;
track.Initialise(&gas, verbose);
track.EnablePlotting(&driftView);

// total number of electrons
int nsum = 0;
// number of attached electrons
int asum = 0;
// energy loss
double esum = 0;
// Initial position and direction 
double x0 = 0., y0 = 0., z0 = 0., t0 = 0.;
double dx0 = 1., dy0 = 0., dz0 = 0.; 
track.NewTrack(x0, y0, z0, t0, dx0, dy0, dz0);

for (const auto& cluster : track.GetClusters()) {
  esum += cluster.energy;
  for (const auto& electron : cluster.electrons) {
    aval.AvalancheElectron(electron.x, electron.y, electron.z, electron.t, 0.1, 0, 0, 0);
    // Loop over the electron avalanches.
    for (int j = 0; j < aval.GetNumberOfElectronEndpoints(); ++j){
      nsum ++;
      double xe1, ye1, ze1, te1, e1;
      double xe2, ye2, ze2, te2, e2;
      int status;
      aval.GetElectronEndpoint(j, xe1, ye1, ze1, te1, e1, xe2, ye2, ze2, te2, e2, status);
      aval.DriftElectron(xe1, ye1, ze1, te1, .1, 0, 0, 0);
      std::cout << "status: " << status << "\n";
      drift.DriftIon(xe1, ye1, ye1, te1);
      if (status == -7) {
        asum ++;
        drift.DriftNegativeIon(xe2, ye2, ye2, te2);
      }
    }
  }
}
double attachRatio = (asum / (double)nsum);


// Print the results
std::cout << "\n" << "Gas: "<< "\n";
for (int i = 0; i < gas.GetNumberOfComponents(); ++i) {
  std::string name;
  double f;
  gas.GetComponent(i, name, f);
  
  std::cout << "  " << name << " : " << f*100 << "\n";
}
std::cout << "  " << gas.GetTemperature()-273.15 << " [°C]" << "\n";
std::cout << "  " << gas.GetPressure() << " [Torr]" << "\n\n";

std::cout << "Gap: " << gap << " [cm]" << "\n";
std::cout << "Voltage: " << voltage << " [V]" << "\n";
std::cout << "Field: " << field << " [V/cm]" << "\n\n";

std::cout << "Energy of particle: " << track.GetEnergy() << " eV" << "\n";
std::cout << "Total energy loss: " << esum << " eV" << "\n\n";

std::cout << "Total electrons: " << nsum << "\n";
std::cout << "Attached electrons: " << asum << "\n";
std::cout << "Ratio of attached electrons: " << attachRatio << "\n";
std::cout << "Free electron fraction: " << 1 - attachRatio << "\n";



if (plotDrift) {
  constexpr bool twod = true;
  constexpr bool drawaxis = true;
  driftView.Plot(twod, drawaxis);
}

Hi,

Adding in the loop @hschindl.

Best,
D

Hi,

The plot function draws straight line segments from the starting point of the track to the first cluster, then from the first cluster to the second cluster and so on. So the plot stops at the last cluster inside the chamber. It doesn’t necessarily mean that the track stops there, but rather that the next cluster would be outside the geometry. I can try to modify the function such that it also plots a line connecting the last cluster inside the geometry to the exit point of the particle.

There is a bit of a flaw in your loop. Every time you call AvalancheElectron or DriftElectron the drift lines stored in AvalancheMicroscopic are reset. So, when you call aval.DriftElectron(...) inside the loop over the electron endpoints, the electron endpoints are reset.
In fact, there is no need to call DriftElectron at all inside the body of the for loop. All electron trajectories that are part of the avalanche are already simulated within AvalancheElectron.

Not right now, but I will implement a method to do that. Should be straightforward, will let you know when it’s done.

Hello @hschindl , thank you so much for all these clarifications. I removed aval.DriftElectron(...) and now there is no more problem.

PS: you should now be able to set the colours of positive and negative ions separately:

1 Like