Simulate clusters of a track

I want to get the clusters information of a track by using trackHeed, but it seems that Heed does not take the electric and magnetic field in the sensor into account for calculating the charged-particle trajectory.
I tried to use the EnableElectricField() &EnableMagneticField() function ,but it didn’t work.

If I want to simulate the relationship between the number of clusters and the electric and magnetic fields, what should I do

Hi @zhong; I am sure @hschindl can help you.


OK,Thanks a lot! :wink:

Can you show your code?

In principle, Heed should take the electric/magnetic field into account when computing the trajectory but it might need a bit of fiddling with the stepping parameters.

Can you provide a minimal working example?

Here’s my working code.

#include <iostream>
#include <fstream>
#include <string>
#include <sstream>
#include <iomanip>
#include <cmath>
#include <vector>

#include <TCanvas.h>
#include <TMath.h>
#include <TROOT.h>
#include <TTree.h>
#include <TApplication.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TBuffer3D.h>
#include <TBuffer3DTypes.h>
#include <TVirtualViewer3D.h>
#include <TFile.h>
#include <TApplication.h>

#include "Garfield/MediumMagboltz.hh"
#include "Garfield/ComponentAnsys123.hh"
#include "Garfield/Sensor.hh"
#include "Garfield/DriftLineRKF.hh"
#include "Garfield/AvalancheMicroscopic.hh"
#include "Garfield/AvalancheMC.hh"
#include "Garfield/ViewGeometry.hh"
#include "Garfield/Plotting.hh"
#include "Garfield/ViewField.hh"
#include "Garfield/ViewFEMesh.hh"
#include "Garfield/ViewSignal.hh"
#include "Garfield/TrackHeed.hh"
#include "Garfield/TrackSrim.hh"

#include "Garfield/SolidBox.hh"
#include "Garfield/GeometrySimple.hh"
#include "Garfield/ComponentAnalyticField.hh"
#include <time.h>

using namespace Garfield;
using namespace std;
using namespace ROOT;

int main() {
    bool debugging = true;
    double trackMom = 100.;
    int nTrk = 100;

    //save cluster information
    double xc = 0., yc = 0., zc = 0., tc = 0.;
    int nCluster;
    double trkion;
    //int Clus_posi;

    TFile *file = new TFile("ar_90_co2_10_0T_1GeV.root","RECREATE");
    TTree* c_position = new TTree("c","c Position");
    //set gas information
    MediumMagboltz* gas = new MediumMagboltz();

    //gas->EnablePenningTransfer(0.44, 0.0, "He");

    // geo
    double zlength = 50.;
    SolidBox* box = new SolidBox(0., 0., 0., 20., 20, zlength); // center and half length
    GeometrySimple* geo = new GeometrySimple();
    geo->AddSolid(box, gas);

    ComponentAnalyticField* comp = new ComponentAnalyticField();
    if(debugging) comp->EnableDebugging();
    else comp->DisableDebugging();
    comp->AddPlaneY(-20, 0., "pad_plane");
    comp->AddPlaneY(20, 10000., "HV");
    comp->SetMagneticField(0, 1, 0);

    Sensor* sensor = new Sensor();
    if(debugging) sensor->EnableDebugging();

    //setup a track 

    double x0 = 0, y0 =0 , z0 = -1*zlength, t0 = 0.; 
    //double dx0 = 0., dy0 = 1., dz0 = dy0/tan(theta);
    double dx0 = 0., dy0 = 0., dz0 = 1.;
    //double dy0 = sin(3.1415926*3/180) , dx0 = 0, dz0 = cos(3.1415926*3/180);
    TrackHeed* track = new TrackHeed();

    // track->EnableElectricField();
    // track->EnableMagneticField();
    // double maxrange = 0., rforstraight = 0., stepstraight = 0., stepcurved = 0.;
    // track->GetSteppingLimits(maxrange, rforstraight, stepstraight, stepcurved);
    // stepcurved = 0.02;
    // track->SetSteppingLimits(maxrange, rforstraight, stepstraight, stepcurved);

    for(int iTrk=0; iTrk<nTrk; iTrk++){
        track->NewTrack(x0, y0, z0, t0, dx0, dy0, dz0);
        // Cluster coordinates
        //double xc = 0., yc = 0., zc = 0., tc = 0.;
        // Number of electrons produced in a collision
        int nc = 0;
        // Energy loss in a collision
        double ec = 0.;
        // Dummy variable (not used at present)
        double extra = 0.;
        //reset cluster
        nCluster = 0;
        //loop clusters 
        while(track->GetCluster(xc, yc, zc, tc, nc, ec, extra)){
            if((zc>= -50.) && (zc< 50.)){
            //Clus_posi = zc;
	        //cout << "cluster position : "<< zc << " deposited energy : "<< ec << endl;
        //cout<<"Cluster number: " << nCluster/100. <<endl;
        //trkion = nCluster/100.0;
    delete track;
    cout<<"end of the program"<<endl;

yes, I have shown it, looking forward to your help

At the momentum (100 GeV/c) and magnetic field (1 T) you are using, you wouldn’t expect a significant curvature of the track over the length of your detector (you can use the formula
p [GeV/c] = 0.3 B [T] r [m] to calculate the curvature radius).
If you lower the momentum you should be able to see a curved track (see the attached example).

test.C (1.3 KB)

Alright, Thanks for your help.

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.