当前位置：文档之家› GEANT4中AIDA的调用

GEANT4中AIDA的调用

Geant4 学习心得—AIDA的使用I

AIDA(Abstract Interfaces for Data Analysis)主要是提供Hbook，ROOT，或者XML格式输出。在学习的过程中我们看到好几种不同的调用方式，在这里我们逐个总结出来。

extended\electromagnetic\TestEM4例子的AIDA调用方式：

（1）在https://www.doczj.com/doc/6d14515248.html,中定义hbook 一维谱（或者root，XML格式的文件）

#ifdef G4ANALYSIS_USE

#include "AIDA/AIDA.h" //需要的头文件

#endif

RunAction::RunAction():af(0), tree(0)

{

histo[0] = 0;

#ifdef G4ANALYSIS_USE

// Creating the analysis factory

af = AIDA_createAnalysisFactory();

if (af) {

// Creating the tree factory

AIDA::ITreeFactory* tf = af->createTreeFactory();

// Creating a tree mapped to an hbook file.

G4bool readOnly = false;

G4bool createNew = true;

G4String options = "--noErrors uncompress";

tree = tf->create("TestEm4.hbook","hbook",readOnly,createNew, options);

//tree = tf->create("TestEm4.root","root",readOnly,createNew, options);

//tree = tf->create("TestEm4.XML" ,"XML" ,readOnly,createNew, options);

delete tf;

if (tree) {

// Creating a histogram factory

AIDA::IHistogramFactory* hf = af->createHistogramFactory(*tree);

// Creating the histogram

histo[0]=hf->createHistogram1D // 头文件中定义了histo[1]

("1","total energy deposit in HpGe Detector(keV)",3000,1.,3000);

// energy ragin 1-5000keV

delete hf;

}

#endif

}

RunAction::~RunAction()

#ifdef G4ANALYSIS_USE

tree->commit(); // Writing the histograms to the file

tree->close(); // and closing the tree (and the file)

delete tree;

delete af;

#endif

}

（2）在https://www.doczj.com/doc/6d14515248.html,中累积谱

#ifdef G4ANALYSIS_USE

#include "AIDA/IHistogram1D.h"

#endif

void EventAction::BeginOfEventAction( const G4Event* evt)

{

G4int evtNb = evt->GetEventID();

//additional initializations

TotalEnergyDeposit = 0.; // 初始化能量变量

}

void EventAction::EndOfEventAction( const G4Event* evt)

{

#ifdef G4ANALYSIS_USE

Run->GetHisto(0)->fill(TotalEnergyDeposit/MeV);

#endif

}

总结：这是最简单的定义一维Histogram的办法。

Geant4 学习心得—AIDA的使用II

AIDA(Abstract Interfaces for Data Analysis)主要是提供Hbook，ROOT，或者XML格式输出。在学习的过程中有几种不同的调用方式。

extended\hadronic\Hadr01例子的AIDA调用方式：

（1）在https://www.doczj.com/doc/6d14515248.html,中调用HistoManager管理器

void RunAction::EndOfRunAction(const G4Run*)

{

HistoManager::GetPointer()->EndOfRun(); //获得Histo的指针}

（2）在https://www.doczj.com/doc/6d14515248.html,中定义HisyoManger的子函数等HistoManager* HistoManager::fManager = 0;

HistoManager* HistoManager::GetPointer()

{

if(!fManager) {

static HistoManager manager;

fManager = &manager;

}

return fManager;

}

（3）在https://www.doczj.com/doc/6d14515248.html,中定义Hbook文件等

#include "Histo.hh"

#ifdef G4ANALYSIS_USE

#include

#include "HistoMessenger.hh"

#endif

Histo::Histo(G4int ver)

{

verbose = ver;

histName = "histo";

histType = "hbook";

option = "--noErrors uncompress";

nHisto = 0;

defaultAct = 1;

tupleName = "tuple.hbook";

tupleId = "100";

tupleList = "";

ntup = 0;

messenger = 0;

#ifdef G4ANALYSIS_USE

messenger = new HistoMessenger(this); // messange 指针tree = 0;

af = 0;

#endif

}

Histo::~Histo()

{

#ifdef G4ANALYSIS_USE

delete messenger;

delete af;

#endif

}

void Histo::book()

{

#ifdef G4ANALYSIS_USE

// Creating the analysis factory

if(!af) af = AIDA_createAnalysisFactory();

// Creating the tree factory

AIDA::ITreeFactory* tf = af->createTreeFactory();

// Creating a tree mapped to a new hbook file.

G4String name = histDir + histName + histExt + "." + histType;

tree = tf->create(name, histType, false, true, option);

delete tf;

// Creating a histogram factory, whose histograms will be handled by the tree AIDA::IHistogramFactory* hf = af->createHistogramFactory( *tree );

// Creating an 1-dimensional histograms in the root directory of the tree

for(G4int i=0; i

if(active) {

G4String idd;

if(histType == "root") idd = "h" + ids;

else idd = ids;

histo = hf->createHistogram1D(idd, tittles, bins, xmin, xmax);

} else {

histo = 0;

}

delete hf;

// Creating a tuple factory, whose tuples will be handled by the tree

if(tupleList != "") {

AIDA::ITupleFactory* tpf = af->createTupleFactory( *tree );

ntup = tpf->create(tupleId, tupleName, tupleList);

delete tpf;

}

#endif

}

void Histo::save()

{

#ifdef G4ANALYSIS_USE

// Write histogram file

tree->commit();

tree->close();

delete tree;

tree = 0;

#endif

}

void Histo::reset()

{

#ifdef G4ANALYSIS_USE

delete tree;

tree = 0;

#endif

}

void Histo::setFileType(const G4String& nam)

{

if(nam == "hbook" || nam == "root" || nam == "aida") histType = nam; else if(nam == "XML" || nam == "xml") histType = "aida";

}

void Histo::add1D(const G4String& id, const G4String& name, G4int nb,

G4double x1, G4double x2, G4double u)

{

if(nHisto > 0) {

for(G4int i=0; i

if(ids == id) return;

}

nHisto++;

x1 /= u;

x2 /= u;

active.push_back(defaultAct);

bins.push_back(nb);

xmin.push_back(x1);

xmax.push_back(x2);

unit.push_back(u);

ids.push_back(id);

tittles.push_back(name);

histo.push_back(0);

}

void Histo::setHisto1D(G4int i, G4int nb, G4double x1, G4double x2, G4double u) {

if(i>=0 && i

bins = nb;

xmin = x1;

xmax = x2;

unit = u;

}

void Histo::fill(G4int i, G4double x, G4double w)

{

#ifdef G4ANALYSIS_USE

if(i>=0 && i

histo->fill(x/unit, w);

}

#endif

}

void Histo::scale(G4int i, G4double x)

{

#ifdef G4ANALYSIS_USE

if(i>=0 && i

histo->scale(x);

}

#endif

}

void Histo::addTuple(const G4String& w1, const G4String& w2, const G4String& w3) {

tupleId = w1;

tupleName = w2;

tupleList = w3;

}

void Histo::fillTuple(const G4String& parname, G4double x)

{

#ifdef G4ANALYSIS_USE

if(ntup) ntup->fill(ntup->findColumn(parname), (float)x);

#endif

}

void Histo::addRow()

{

#ifdef G4ANALYSIS_USE

if(ntup) ntup->addRow();

#endif

}

void Histo::print(G4int i)

{

#ifdef G4ANALYSIS_USE

if(i>=0 && i

G4double step = (xmax - xmin)/G4double( bins);

G4double x = xmin - step*0.5;

G4double y, maxX=0, maxY=0;

G4int maxJ=0;

for(G4int j=0; j

x += step;

y = histo->binHeight(j);

if(maxY < y) {maxY = y; maxX = x; maxJ = j;}

}

#endif

}

（4）在RUNACTION中开始和结束累谱。

void RunAction::BeginOfRunAction(const G4Run* aRun)

{

G4int id = aRun->GetRunID();

(HistoManager::GetPointer())->BeginOfRun();

void RunAction::EndOfRunAction(const G4Run*)

{

HistoManager::GetPointer()->EndOfRun(); //开始累谱

}

（5）在灵敏探测器的定义中

TargetSD::TargetSD(const G4String& name)

:G4VSensitiveDetector(name)

{

theHisto = HistoManager::GetPointer();

}

G4bool TargetSD::ProcessHits(G4Step* aStep, G4TouchableHistory*)

{

theHisto->AddTargetStep(aStep);

return true;

}

CheckVolumeSD::CheckVolumeSD(const G4String& name)

:G4VSensitiveDetector(name)

{

theHisto = HistoManager::GetPointer();

}

G4bool CheckVolumeSD::ProcessHits(G4Step* aStep, G4TouchableHistory*)

{

const G4Track* track = aStep->GetTrack();

if(track->GetTrackID() > 1) theHisto->AddLeakingParticle(track);

return true;

Geant4 学习心得—AIDA的使用III

AIDA(Abstract Interfaces for Data Analysis)主要是提供Hbook，ROOT，或者XML格式输出。在学习的过程中有几种不同的调用方式。

extended\biasing\B02例子的AIDA调用方式：

（1）在https://www.doczj.com/doc/6d14515248.html,中调用AIDA

// AIDA stuff

#include "AIDA/IAnalysisFactory.h"

#include "AIDA/ITree.h"

#include "AIDA/ITreeFactory.h"

#include "AIDA/IHistogram1D.h"

#include "AIDA/IHistogramFactory.h"

int main(int , char **)

{

// create a histogram for a special scorer for the last cell

AIDA::IAnalysisFactory *af = AIDA_createAnalysisFactory(); AIDA::ITreeFactory *tf = af->createTreeFactory();

AIDA::ITree *tree = tf->create("b02.hbook", "hbook",false,true); AIDA::IHistogramFactory *hf = af->createHistogramFactory( *tree ); AIDA::IHistogram1D *h =

hf->createHistogram1D("10","w*sl vs. e", 30, 0., 20*MeV);

tree->commit();

tree->close();

delete af;

delete tf;

delete tree;

delete runManager;

return 0;

}

（2）在Score中调用Histo，在hh头文件中定义了fHisto

B02CellScorer::B02CellScorer(AIDA::IHistogram1D *h):fHisto(h){}

void B02CellScorer::ScoreAnExitingStep(const G4Step &aStep,

const G4GeometryCell &pre_gCell)

{

fG4CellScorer.ScoreAnExitingStep(aStep, pre_gCell);

FillHisto(aStep);

}

void B02CellScorer::FillHisto(const G4Step &aStep)

{

G4StepPoint *p = 0;

p = aStep.GetPreStepPoint();

G4double sl = aStep.GetStepLength();

G4double slw = sl * p->GetWeight();

fHisto->fill(p->GetKineticEnergy(), slw);

}

Geant4 学习心得—AIDA的使用IV

AIDA(Abstract Interfaces for Data Analysis)主要是提供Hbook，ROOT，或者XML格式输出。在学习的过程中有几种不同的调用方式。

extended\Analysis\AnaEx01例子的AIDA调用方式：采用https://www.doczj.com/doc/6d14515248.html,

（1）在https://www.doczj.com/doc/6d14515248.html,中定义所有的Histogram，Tuple 以及Plotter等

#ifdef G4ANALYSIS_USE

#include "G4ios.hh"

#include "G4SDManager.hh"

#include "G4Run.hh"

#include "G4Event.hh"

#include "G4HCofThisEvent.hh"

#include

#include "AnaEx01CalorHit.hh"

#include "AnaEx01AnalysisManager.hh"

AnaEx01AnalysisManager::AnaEx01AnalysisManager(AIDA::IAnalysisFactory* aAIDA) :fCalorimeterCollID(-1),fAIDA(aAIDA),fTree(0),fEAbs(0),fLAbs(0),fEGap(0),

fLGap(0),fTuple(0)

{

// Could fail if no AIDA implementation found :

if(!fAIDA) {

G4cout << "AIDA analysis factory not found." << G4endl;

return;

}

AIDA::ITreeFactory* treeFactory = fAIDA->createTreeFactory();

if(!treeFactory) return;

// Create a tree-like container to handle histograms.

// This tree is associated to a AnaEx01. file.

std::string h_name_EAbs("EAbs");

std::string h_name_LAbs("LAbs");

std::string h_name_EGap("EGap");

std::string h_name_LGap("LGap");

std::string t_name("AnaEx01");

// File format :

std::string format("xml");

//std::string format("hbook");

//std::string format("root");

std::string file("AnaEx01");

std::string ext;

ext = "."+format;

std::string opts;

if(format=="hbook") {

opts = "compress=no";

h_name_EAbs = "1";

h_name_LAbs = "2";

h_name_EGap = "3";

h_name_LGap = "4";

t_name = "101";

} else if(format=="root") {

opts = "compress=yes";

} else if(format=="xml") {

ext = ".aida";

opts = "compress=no";

} else {

G4cout << "storage format \"" << format << "\""

<< " not handled in this example."

<< G4endl;

return;

}

file += ext;

fTree = treeFactory->create(file,format,false,true,opts);

// Factories are not "managed" by an AIDA analysis system. // They must be deleted by the AIDA user code.

delete treeFactory;

if(!fTree) {

G4cout << "can't create tree associated to file \"" << file << "\"."

<< G4endl;

return;

}

fTree->mkdir("histograms");

fTree->cd("histograms");

// Create an histo factory that will create histo in the tree :

AIDA::IHistogramFactory* histoFactory =

fAIDA->createHistogramFactory(*fTree);

if(histoFactory) {

fEAbs = histoFactory->createHistogram1D(h_name_EAbs,"EAbs",100,0,100);

if(!fEAbs) G4cout << "can't create histo EAbs." << G4endl;

fLAbs = histoFactory->createHistogram1D(h_name_LAbs,"LAbs",100,0,100);

if(!fLAbs) G4cout << "can't create histo LAbs." << G4endl;

fEGap = histoFactory->createHistogram1D(h_name_EGap,"EGap",100,0,10);

if(!fEGap) G4cout << "can't create histo EGap." << G4endl;

fLGap = histoFactory->createHistogram1D(h_name_LGap,"LGap",100,0,100);

if(!fLGap) G4cout << "can't create histo LGap." << G4endl;

delete histoFactory;

}

fTree->cd("..");

fTree->mkdir("tuples");

fTree->cd("tuples");

// Get a tuple factory :

AIDA::ITupleFactory* tupleFactory = fAIDA->createTupleFactory(*fTree);

if(tupleFactory) {

// Create a tuple :

fTuple = tupleFactory->create(t_name,"AnaEx01",

"double EAbs,double LAbs,double EGap,double LGap");

if(!fTuple) G4cout << "can't create tuple." << G4endl;

delete tupleFactory;

}

fTree->cd("..");

}

AnaEx01AnalysisManager::~AnaEx01AnalysisManager() {

}

void AnaEx01AnalysisManager::BeginOfRun(const G4Run* aRun){

G4cout << "### Run " << aRun->GetRunID() << " start." << G4endl;

}

void AnaEx01AnalysisManager::EndOfRun(const G4Run*){

if(fTree) fTree->commit();

if(fEAbs) {

G4cout << "Histo : EAbs : mean " << fEAbs->mean() << " rms : " << fEAbs->rms() << G4endl;

G4cout << "Histo : LAbs : mean " << fLAbs->mean() << " rms : " << fLAbs->rms() << G4endl;

G4cout << "Histo : EGap : mean " << fEGap->mean() << " rms : " << fEGap->rms() << G4endl;

G4cout << "Histo : LGap : mean " << fLGap->mean() << " rms : " << fLGap->rms() << G4endl;

}

void AnaEx01AnalysisManager::BeginOfEvent(const G4Event*){

if(fCalorimeterCollID==-1) {

G4SDManager* SDman = G4SDManager::GetSDMpointer();

fCalorimeterCollID = SDman->GetCollectionID("CalCollection");

}

void AnaEx01AnalysisManager::EndOfEvent(const G4Event* aEvent){

if(!fEAbs) return; // No histo booked !

if(!fTuple) return; // No tuple booked !

//G4int evtNb = aEvent->GetEventID();

G4HCofThisEvent* HCE = aEvent->GetHCofThisEvent();

AnaEx01CalorHitsCollection* CHC =

HCE ? (AnaEx01CalorHitsCollection*)(HCE->GetHC(fCalorimeterCollID)) : 0;

if(CHC) {

G4int n_hit = CHC->entries();

for (G4int i=0;i

G4double EAbs = (*CHC)->GetEdepAbs();

G4double LAbs = (*CHC)->GetTrakAbs();

G4double EGap = (*CHC)->GetEdepGap();

G4double LGap = (*CHC)->GetTrakGap();

fEAbs->fill(EAbs);

fLAbs->fill(LAbs);

fEGap->fill(EGap);

fLGap->fill(LGap);

fTuple->fill(0,EAbs);

fTuple->fill(1,LAbs);

fTuple->fill(2,EGap);

fTuple->fill(3,LGap);

fTuple->addRow();

}

void AnaEx01AnalysisManager::Step(const G4Step*){} #endif

（2）头文件定义为：

#ifndef AnaEx01AnalysisManager_h

#define AnaEx01AnalysisManager_h 1

#ifdef G4ANALYSIS_USE

class G4Run;

class G4Event;

class G4Step;

namespace AIDA {

class IAnalysisFactory;

class ITree;

class IHistogram1D;

class ITuple;

}

class AnaEx01AnalysisManager {

public:

AnaEx01AnalysisManager(AIDA::IAnalysisFactory*); virtual ~AnaEx01AnalysisManager();

public:

virtual void BeginOfRun(const G4Run*);

virtual void EndOfRun(const G4Run*);

virtual void BeginOfEvent(const G4Event*);

virtual void EndOfEvent(const G4Event*);

virtual void Step(const G4Step*);

private:

int fCalorimeterCollID;

AIDA::IAnalysisFactory* fAIDA;

AIDA::ITree* fTree;

AIDA::IHistogram1D* fEAbs;

AIDA::IHistogram1D* fLAbs;

AIDA::IHistogram1D* fEGap;

AIDA::IHistogram1D* fLGap;

AIDA::ITuple* fTuple;

};

#else

class AnaEx01AnalysisManager;

#endif

（3）在RunAction中调用，填谱

#ifdef G4ANALYSIS_USE

#include "AnaEx01AnalysisManager.hh"

#endif

#include "AnaEx01RunAction.hh"

AnaEx01RunAction::AnaEx01RunAction(

AnaEx01AnalysisManager* aAnalysisManager

):fAnalysisManager(aAnalysisManager){}

AnaEx01RunAction::~AnaEx01RunAction(){}

void AnaEx01RunAction::BeginOfRunAction(const G4Run* aRun) { #ifdef G4ANALYSIS_USE

if(fAnalysisManager) fAnalysisManager->BeginOfRun(aRun);

#endif

}

void AnaEx01RunAction::EndOfRunAction(const G4Run* aRun){ #ifdef G4ANALYSIS_USE

if(fAnalysisManager) fAnalysisManager->EndOfRun(aRun);

#endif

}

Geant4 学习心得—如何转动实体

我们通过A01例子来学习如何设置实体的转动

(1) 在DetectorConstructtion.hh中定义RotationMatrix：

#include "G4RotationMatrix.hh"

class A01DetectorConstruction : public G4VUserDetectorConstruction

{

private:

G4RotationMatrix* armRotation;

}

（2）在https://www.doczj.com/doc/6d14515248.html,中定义转动角度，在G4VPlacement中的第一项设置转动。armRotation = new G4RotationMatrix();

armRotation->rotateY(armAngle);

// second arm

secondArmPhys

= new G4PVPlacement(armRotation,G4ThreeVector(x,0.,z),secondArmLogical,

"secondArmPhys",worldLogical,0,0);

（3）A01中https://www.doczj.com/doc/6d14515248.html,的设置磁场转动

// placement of Tube

G4RotationMatrix* fieldRot = new G4RotationMatrix();

fieldRot->rotateX(90.*deg);

new G4PVPlacement(fieldRot,G4ThreeVector(),magneticLogical,

"magneticPhysical",worldLogical,0,0);

(4) 在ExN04的https://www.doczj.com/doc/6d14515248.html,中设置的转动

G4RotationMatrix rm;

rm.rotateZ(phi);

muoncounter_phys = new G4PVPlacement(G4Transform3D(rm,G4ThreeVector(x,y,z)),

muoncounter_log, "muoncounter_P",

experimentalHall_log,false,i);