#include #include #include #include "math.h" #include #include #include #include #include #include "Garfield/MediumConductor.hh" #include "Garfield/ViewGeometry.hh" #include "Garfield/SolidBox.hh" #include "Garfield/ComponentAnsys123.hh" #include "Garfield/ComponentComsol.hh" #include "Garfield/ViewField.hh" #include "Garfield/ViewFEMesh.hh" #include "Garfield/MediumMagboltz.hh" #include "Garfield/Sensor.hh" #include "Garfield/AvalancheMicroscopic.hh" #include "Garfield/AvalancheMC.hh" #include "Garfield/Random.hh" #include "Garfield/TrackSrim.hh" #include "Garfield/GeometrySimple.hh" #include "Garfield/ViewSignal.hh" #include "Garfield/ViewCell.hh" #include "Garfield/Plotting.hh" #include #include #include #include #include #include #include #include #include "Garfield/ComponentAnsys123.hh" #include "Garfield/ViewField.hh" #include "Garfield/ViewFEMesh.hh" #include "Garfield/MediumMagboltz.hh" #include "Garfield/Sensor.hh" #include "Garfield/AvalancheMicroscopic.hh" #include "Garfield/AvalancheMC.hh" #include "Garfield/Random.hh" using namespace Garfield; using namespace std; void ModifyLineData(string fileName, int lineNum, char* lineData) { cout<<"ModifyLineData"<SetComponent(&fm); fieldViewXZ->SetPlane(0, -0.1, 0, 0, 0.05, 0); fieldViewXZ->SetArea(-Dim_x, Dim_z_down, Dim_x, Dim_z_up); fieldViewXZ->SetSensor(&sensor); fieldViewXZ->SetNumberOfSamples2d(500,500); //fieldViewXZ->SetElectricFieldRange(0,2000); fieldViewXZ->SetCanvas(cFieldXZ); fieldViewXZ->Plot("e","colz"); cFieldXY = new TCanvas("cFieldXY", "", 900, 800); fieldViewXY->SetComponent(&fm); fieldViewXY->SetPlane(0, 0, -0.1, 0, 0, 0.05); fieldViewXY->SetArea(-Dim_x, -Dim_y, Dim_x, Dim_y); fieldViewXY->SetSensor(&sensor); fieldViewXY->SetNumberOfSamples2d(500,500); //fieldViewXY->SetElectricFieldRange(0,2000); fieldViewXY->SetCanvas(cFieldXY); fieldViewXY->Plot("e","colz"); //fieldViewXY->PlotWeightingField(label_Anode1,"p","colz"); } // // ---- Track ---- // TrackSrim *tr = new TrackSrim(); tr->SetSensor(&sensor); tr->ReadFile("./gasfile/alpha.txt"); // 导入气体文件 tr->SetKineticEnergy(5.486e6); // 粒子能量5.486MeV tr->SetWorkFunction(26580); // 平均电离能 , 为了计算速度加快,在此放大了1000倍 tr->SetFanoFactor(0.3); // 法诺因子 constexpr double za = 0.9 * (18. / 40.) + 0.1 * (10. / 16.); tr->SetAtomicMassNumbers(10. / za, 10.); // 气体靶的A/Z值 tr->SetClustersMaximum(100); // 一个团簇多少电子 tr->Print(); // // ---- Micro 电子漂移轨迹 ---- // //AvalancheMicroscopic *aval = new AvalancheMicroscopic(); AvalancheMC*aval = new AvalancheMC(); //aval->EnableWeightingFieldIntegration(); aval->SetSensor(&sensor); //aval->EnableAttachmentMarkers(); aval->EnableSignalCalculation(true); aval->UseWeightingPotential(true); //aval->SetTimeWindow(tMin, tMax); //aval->SetElectronTransportCut(0.); aval->DisableAvalancheSizeLimit(); //aval->EnableDriftLines(); // // ---- 几何、track、电子轨迹 可视化 ---- // /* //ViewCell cellView; TCanvas canvas("cCanvas", "", 600, 600); canvas.SetLeftMargin(0.16); */ ViewDrift *driftView = new ViewDrift(); if (plotDrift) { cDrift = new TCanvas("cDrift", "", 600, 600); driftView->SetClusterMarkerSize(0.1); driftView->SetCollisionMarkerSize(0.1); driftView->SetColourElectrons(kRed); driftView->SetColourTracks(kGreen + 3); tr->EnablePlotting(driftView); aval->EnablePlotting(driftView); //driftView->SetCanvas(&canvas); } TRandom rrand; double rand; const unsigned int nTracks = 30; // --- 释放的Alpha 粒子个数 --- // TLatex label; int count = 0; for (unsigned int j = 0; j < nTracks; ++j) { //sensor.ClearSignal(); cout << " ------------------------ New Track ------------------------ "<NewTrack(x0, y0, z0, t0, dx, dy, dz); // loop clusters while(tr->GetCluster(xc, yc, zc, tc, nc, ec, ekin)) //获取当前团簇信息 { ee1 = 0.1; //当前cluster的平均电子初始能量 cout << " ===================================================== "<GetAvalancheSize(ne, ni); //雪崩电子数,(不雪崩) //driftView->Clear(); aval->DriftElectron(xc, yc, zc, tc); //开始计算电子,ee1,初电子各向同性 const int np = aval->GetNumberOfElectronEndpoints(); //倍增电子个数,1 cout << " Multiplicity: " << np << endl; for(int j = 0; j < np; j++) { aval->GetElectronEndpoint(j, xe1, ye1, ze1, te1, xe2, ye2, ze2, te2, status); //倍增电子最终状态 cout << " Electron start position: " << xe1 <<","<< ye1 <<","<< ze1 <<"\t"<< " Electron stop position: " << xe2 <<","<< ye2 <<","<< ze2 <Fill(); Total_counts++; } } } } // // ---- 几何、track、电子轨迹 可视化 ---- // if (plotDrift) { constexpr bool plotMesh = false; if (plotMesh) { ViewFEMesh* meshView = new ViewFEMesh(); meshView->SetArea(-Dim_x/20, -Dim_y/20, Dim_z_down, Dim_x/20, Dim_y/20, Dim_z_up); meshView->SetCanvas(cDrift); meshView->SetComponent(&fm); // x-z projection. meshView->SetPlane(0, -0.1, 0, 0, 0.05, 0); meshView->SetFillMesh(true); driftView->SetColourElectrons(kRed); driftView->SetColourTracks(kGreen + 3); // Set the color. meshView->SetViewDrift(driftView); meshView->Plot(); } else { driftView->SetPlane(0, -0.1, 0, 0, 0, 0); driftView->SetArea(-Dim_x/4., Dim_z_down/4., Dim_x/4., Dim_z_up/4.); driftView->SetCanvas(cDrift); constexpr bool twod = true; driftView->Plot(twod, twod); } } ofstream signal_pAnode1; signal_pAnode1.open( "SignalpAnode1.txt"); ofstream signal_pAnode2; signal_pAnode2.open( "SignalpAnode2.txt"); ofstream signal_pCathode; signal_pCathode.open( "SignalpCathode.txt"); double Q_Anode1=0; double Q_Anode2=0; double Q_Cathode=0; for(int i = 0; i < nTimeBins; i++) { Q_Anode1 = (sensor.GetSignal(label_pAnode1, i)); signal_pAnode1 << i << "\t" << Q_Anode1 << endl; Q_Anode2 = (sensor.GetSignal(label_pAnode2, i)); signal_pAnode2 << i << "\t" << Q_Anode2 << endl; Q_Cathode = (sensor.GetSignal(label_pCathode, i)); signal_pCathode << i << "\t" << Q_Cathode << endl; } // // ---- 输出电荷信号可视化 ---- // ViewSignal Signalview_pAnode1; // pAnode1 信号可视化 ViewSignal Signalview_pAnode2; // pAnode2 信号可视化 ViewSignal Signalview_pCathode; if (plotSignal) { ChargeSignal_pAnode1 = new TCanvas("Anode1 Signal", "", 600, 600); Signalview_pAnode1.SetSensor(&sensor); sensor.IntegrateSignal(label_pAnode1); Signalview_pAnode1.SetCanvas(ChargeSignal_pAnode1); Signalview_pAnode1.SetLabelY("Induced Charge [fC]"); Signalview_pAnode1.PlotSignal(label_pAnode1,true,false,false,true,false); ChargeSignal_pAnode2 = new TCanvas("Anode2 Signal", "", 600, 600); Signalview_pAnode2.SetSensor(&sensor); sensor.IntegrateSignal(label_pAnode2); Signalview_pAnode2.SetCanvas(ChargeSignal_pAnode2); Signalview_pAnode2.SetLabelY("Induced Charge [fC]"); Signalview_pAnode2.PlotSignal(label_pAnode2,true,false,false,true,false); ChargeSignal_pCathode = new TCanvas("Cathode Signal", "", 600, 600); Signalview_pCathode.SetSensor(&sensor); sensor.IntegrateSignal(label_pCathode); Signalview_pCathode.SetCanvas(ChargeSignal_pCathode); Signalview_pCathode.SetSensor(&sensor); Signalview_pCathode.SetLabelY("Induced Charge [fC]"); Signalview_pCathode.PlotSignal(label_pCathode); } if(savefile) { string rootfile = "/Results.root"; TFile *hist = new TFile("Results.root", "RECREATE"); ChargeSignal_pAnode1->Write("Anode1 Signal"); ChargeSignal_pAnode2->Write("Anode2 Signal"); ChargeSignal_pCathode->Write("Cathode Signal"); cFieldXY->Write("FieldXY"); cFieldXZ->Write("FieldXZ"); cDrift->Write("Drift"); hist->Close(); } cout << "=======================" <