-Compare.C
-Config.C.all
-Config.C.keep
-Digitize.C
-Document.C
-LsHits.C
-MakeMap.C
-RawTest.C
-ReadRaw.C
-ShowDigits.C
-ShowHits.C
-old
-runIt.C
-runflukageo.sh
-rungeant3geo.sh
-va1_train.C
--- /dev/null
+//
+// Script to compare the output of GEANT 3.21 to FLUKA 2.
+//
+void
+Compare()
+{
+ TFile* fluka_file = TFile::Open("fluka/FMD.Hits.root", "READ");
+ TFile* geant_file = TFile::Open("geant321/FMD.Hits.root", "READ");
+ if (!fluka_file || !geant_file) {
+ std::cerr << "Failed to open one or more of the input files"
+ << std::endl;
+ return;
+ }
+
+
+ fluka_file->cd();
+ gDirectory->cd("Event0");
+ TTree* fluka_tree = static_cast<TTree*>(gDirectory->Get("TreeH"));
+
+ geant_file->cd();
+ gDirectory->cd("Event0");
+ TTree* geant_tree = static_cast<TTree*>(gDirectory->Get("TreeH"));
+
+ if (!fluka_tree || !geant_tree) {
+ std::cerr << "Failed to get one or more of the trees"
+ << std::endl;
+ return;
+ }
+
+ gStyle->SetOptStat(0);
+ gStyle->SetOptTitle(0);
+ gStyle->SetLabelFont(132, "XY");
+ gStyle->SetTitleFont(132, "XY");
+ gStyle->SetTextFont(132);
+ gStyle->SetNdivisions(10, "XY");
+
+ TCanvas* c = new TCanvas("c", "c", 600, 600);
+ c->SetBorderMode(0);
+ c->SetBorderSize(0);
+ c->SetFillColor(0);
+ c->SetTopMargin(.05);
+ c->SetRightMargin(.05);
+
+ TH1* fluka_dist = new TH1F("fluka_dist", "FLUKA Energy deposition",
+ 100, 0, 1);
+ fluka_dist->SetLineColor(2);
+ fluka_dist->SetFillColor(2);
+ fluka_dist->SetFillStyle(3001);
+ fluka_dist->GetXaxis()->SetTitle("Energy deposited");
+ fluka_dist->GetYaxis()->SetTitle("Frequency");
+ fluka_tree->Draw("FMD.fEdep>>fluka_dist");
+
+ TH1* geant_dist = new TH1F("geant_dist", "GEANT Energy deposition",
+ 100, 0, 1);
+ geant_dist->SetLineColor(3);
+ geant_dist->SetFillColor(3);
+ geant_dist->SetFillStyle(3002);
+ geant_tree->Draw("FMD.fEdep>>geant_dist", "", "SAME");
+
+ TLegend* l = new TLegend(.3, .8, .95, .95);
+ l->SetFillColor(0);
+ l->SetBorderSize(1);
+ l->SetTextFont(132);
+ l->AddEntry(fluka_dist, Form("%s - %d counts",
+ fluka_dist->GetTitle(),
+ Int_t(fluka_dist->Integral())), "LF");
+ l->AddEntry(geant_dist, Form("%s - %d counts",
+ geant_dist->GetTitle(),
+ Int_t(geant_dist->Integral())), "LF");
+ l->Draw();
+
+ c->Modified();
+ c->cd();
+ c->Print("fluka_vs_geant321.C");
+}
+
+
+
+
--- /dev/null
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * *
+ * Author: The ALICE Off-line Project. *
+ * Contributors are mentioned in the code where appropriate. *
+ * *
+ * Permission to use, copy, modify and distribute this software and its *
+ * documentation strictly for non-commercial purposes is hereby granted *
+ * without fee, provided that the above copyright notice appears in all *
+ * copies and that both the copyright notice and this permission notice *
+ * appear in the supporting documentation. The authors make no claims *
+ * about the suitability of this software for any purpose. It is *
+ * provided "as is" without express or implied warranty. *
+ **************************************************************************/
+
+/* $Id$ */
+
+// Script to document the FMD code
+void
+Document()
+{
+ gEnv->SetValue("Root.Html.SourceDir", "$(ALICE)/FMD");
+ gEnv->SetValue("Root.Html.OutputDir", "$(ALICE)/FMD/html");
+
+ gSystem->MakeDirectory("$(ALICE)/FMD/html");
+
+ THtml* html = new THtml;
+ html->MakeAll(kFALSE, "AliFMD*");
+ html->Convert("$(ALICE)/FMD/Digitize.C", "Digitize",
+ "FMD/html/src");
+ html->Convert("$(ALICE)/FMD/Reconstruct.C", "Reconstruct",
+ "FMD/html/src");
+ html->Convert("$(ALICE)/FMD/Simulate.C", "Simulate",
+ "FMD/html/src");
+ html->Convert("$(ALICE)/FMD/DrawFMD.C", "DrawFMD",
+ "FMD/html/src");
+ html->Convert("$(ALICE)/FMD/ViewFMD.C", "ViewFMD",
+ "FMD/html/src");
+ html->MakeIndex("AliFMD*");
+
+ std::ofstream index("FMD/html/index.html");
+ html->WriteHtmlHeader(index, "ALICE FMD Code - Index page");
+
+ index << "<h1>ALICE FMD Code</h1>\n"
+ << "<ul>\n"
+ << "<li><a href=\"USER_Index.html\">Classes</a></li>\n"
+ << "<li><a href=\"src/Digitize.C.html\">Digitize script</a></li>\n"
+ << "<li><a href=\"src/Reconstruct.C.html\">Reconstruct script</a></li>\n"
+ << "<li><a href=\"src/Simulate.C.html\">Simulate script</a></li>\n"
+ << "<li><a href=\"src/DrawFMD.C.html\">DrawFMD script</a></li>\n"
+ << "<li><a href=\"src/ViewFMD.C.html\">ViewFMD script</a></li>\n"
+ << "</ul>\n"
+ << std::endl;
+ html->WriteHtmlFooter(index, "", "", "", "");
+ index.close();
+}
+
+//
+// EOF
+//
--- /dev/null
+//
+// Script to draw detail of the FMD
+//
+void DrawDetail()
+{
+ // gAlice->Init("FMD/scripts/ConfigInner.C");
+ gAlice->Init("FMD/scripts/ConfigFmdOnly.C");
+ gMC->Gsatt("*", "seen", -1);
+ gMC->Gsatt("alic", "seen", 0);
+ gROOT->LoadMacro("FMD/ViewFMD.C");
+ gInterpreter->ProcessLine("ViewFMD()");
+ // gROOT->LoadMacro("VZERO/ViewVZERO.C");
+ // gInterpreter->ProcessLine("ViewVZERO()");
+ // gROOT->LoadMacro("START/ViewSTART.C");
+ // gInterpreter->ProcessLine("ViewSTART()");
+ // gROOT->LoadMacro("macros/ViewITS.C");
+ // gInterpreter->ProcessLine("ViewITS()");
+ // gROOT->LoadMacro("FMD/scripts/ViewPIPE.C");
+ // gInterpreter->ProcessLine("ViewPIPE()");
+ // gMC->Gsatt("ITSV", "seen", 1);
+ // gMC->Gsatt("0STR", "seen", 1);
+ // gMC->Gsatt("0STL", "seen", 1);
+ // gMC->Gsatt("0SUP", "seen", 1);
+ // gMC->Gsatt("FMD1", "seen", 1);
+ // gMC->Gsatt("FMD2", "seen", 1);
+ // gMC->Gsatt("FOAC", "seen", 1);
+ // gMC->Gsatt("FOVF", "seen", 1);
+ // gMC->Gsatt("FOVB", "seen", 1);
+ // gMC->Gsatt("FIRG", "seen", 1);
+ // gMC->Gsatt("FIRG", "colo", 1);
+ // gMC->Gsatt("F3IH", "seen", -1);
+ // gMC->Gsatt("F3II", "seen", -1);
+ // gMC->Gsatt("F3IJ", "seen", -1);
+ // gMC->Gsatt("F3IK", "seen", -1);
+ gMC->Gsatt("FMD3", "seen", 1);
+ gMC->Gsatt("FMD3", "colo", 3);
+ gMC->Gsatt("FMD3", "lsty", 1);
+ gMC->Gsatt("F3SL", "colo", 2);
+ gMC->Gsatt("F3SN", "colo", 2);
+ gMC->Gsatt("F3SB", "colo", 2);
+ gMC->Gsatt("F3SF", "colo", 2);
+ // gMC->Gsatt("FIAC", "lwid", 2);
+ // gMC->Gsatt("FIAC", "colo", 3);
+ // gMC->Gsatt("FOAC", "lwid", 2);
+ // gMC->Gsatt("FOAC", "colo", 3);
+ gMC->Gdopt("hide", "off");
+ gMC->Gdopt("shad", "off");
+#if 1
+ gMC->Gdopt("hide", "on");
+ // gMC->Gdopt("shad", "on");
+ // gMC->Gsatt("*", "fill", 7);
+#endif
+ gMC->SetClipBox(".");
+ gMC->SetClipBox("*", 0, 10000, -1000, 1000, -1000, 1000);
+ gMC->DefaultRange();
+ gMC->Gdraw("alic", 60, 0, 0, -5, 10, .3, .3);
+
+ gPad->Modified();
+ gPad->cd();
+ // gPad->Print("FMD3_detail.png");
+}
--- /dev/null
+//
+// Script to draw detail of the FMD
+//
+void DrawFMD3()
+{
+ // gAlice->Init("FMD/scripts/ConfigInner.C");
+ gAlice->Init("FMD/scripts/ConfigFmdOnly.C");
+ gMC->Gsatt("*", "seen", -1);
+ gMC->Gsatt("alic", "seen", 0);
+ gROOT->LoadMacro("FMD/ViewFMD.C");
+ gInterpreter->ProcessLine("ViewFMD()");
+ gMC->Gdopt("hide", "on");
+ gMC->Gdopt("shad", "on");
+ gMC->Gsatt("*", "fill", 7);
+ gMC->SetClipBox(".");
+ gMC->SetClipBox("*", 0, 1000, -1000, 1000, -1000, 1000);
+ gMC->DefaultRange();
+ gMC->Gdraw("alic", 60, 0, 0, -3, 10, .25, .25);
+
+ TArrow* a1 = new TArrow(13.5, 16, 15, 18., .03, "<|");
+ a1->SetAngle(30);
+ a1->SetFillColor(1);
+ a1->Draw();
+
+ TLatex* l1 = new TLatex(15, 18, "Honeycomb");
+ l1->SetTextAlign(12);
+ l1->SetTextFont(132);
+ l1->SetTextSize(.04);
+ l1->Draw();
+
+ a1->DrawArrow(13.4, 14., 15, 15, .03, "<|");
+ l1->DrawLatex(15, 15, "Support Leg");
+
+ a1->DrawArrow(10.7, 14.2, 15, 13, .03, "<|");
+ l1->DrawLatex(15, 13, "Print board");
+
+ a1->DrawArrow(9.7, 12.5, 15, 11, .03, "<|");
+ l1->DrawLatex(15, 11, "Silicon sensor");
+
+ a1->DrawArrow(8.3, 12.7, 7, 15, .03, "<|");
+ TLatex* l2 = new TLatex(7, 15, "Support Cone");
+ l2->SetTextSize(.04);
+ l2->SetTextFont(132);
+ l2->SetTextAlign(32);
+ l2->Draw();
+
+ TLatex* l3 = new TLatex(3, 3, "FMD3");
+ l3->SetTextSize(.06);
+ l3->SetTextFont(132);
+ l3->Draw();
+
+ gPad->Modified();
+ gPad->cd();
+ gPad->Print("FMD3_detail.png");
+}
--- /dev/null
+//
+// Script to draw detail of the FMD
+//
+void DrawInner()
+{
+ gAlice->Init("FMD/scripts/ConfigInner.C");
+ gMC->Gsatt("*", "seen", -1);
+ gMC->Gsatt("alic", "seen", 0);
+ gROOT->LoadMacro("FMD/ViewFMD.C");
+ gInterpreter->ProcessLine("ViewFMD()");
+ gROOT->LoadMacro("VZERO/ViewVZERO.C");
+ gInterpreter->ProcessLine("ViewVZERO()");
+ gROOT->LoadMacro("START/ViewSTART.C");
+ gInterpreter->ProcessLine("ViewSTART()");
+ gROOT->LoadMacro("macros/ViewITS.C");
+ gInterpreter->ProcessLine("ViewITS()");
+ // gROOT->LoadMacro("FMD/scripts/ViewPIPE.C");
+ // gInterpreter->ProcessLine("ViewPIPE()");
+ // gMC->Gsatt("ITSV", "seen", 1);
+ gMC->Gsatt("0STR", "seen", 1);
+ gMC->Gsatt("0STL", "seen", 1);
+ gMC->Gsatt("0SUP", "seen", 1);
+ // gMC->Gsatt("FMD1", "seen", 1);
+ // gMC->Gsatt("FMD2", "seen", 1);
+ gMC->Gsatt("FMD3", "seen", 0);
+ gMC->Gdopt("hide", "on");
+ gMC->Gdopt("shad", "on");
+ gMC->Gsatt("*", "fill", 7);
+ gMC->SetClipBox(".");
+ gMC->SetClipBox("*", 0, 1000, -1000, 1000, -1000, 1000);
+ gMC->DefaultRange();
+ gMC->Gdraw("alic", 60, 0, 0, 10, 10, .10, .10);
+ // gMC->Gdhead(1111, "FMD3 detail");
+ // gMC->Gdman(16, 10, "MAN");
+
+ gPad->Modified();
+ gPad->cd();
+ gPad->Print("Inner.png");
+}
--- /dev/null
+//
+// Script to plot the lego information
+//
+#ifndef __CINT__
+# include <TFile.h>
+# include <TCanvas.h>
+# include <TPad.h>
+# include <TH1D.h>
+# include <TString.h>
+# include <TLatex.h>
+# include <TLegend.h>
+# include <THStack.h>
+# include <TStyle.h>
+# include <iostream>
+using namespace std;
+#endif
+
+Float_t
+deg2eta(Float_t ang)
+{
+ if (ang == 180) ang -= .001;
+ if (ang == 0) ang += .001;
+ Float_t theta = ang * TMath::Pi() / 180;
+ Float_t eta = - TMath::Log(TMath::Tan(theta / 2));
+ return eta;
+}
+
+
+TH1*
+getHisto(const Char_t* which, const Char_t* what, TH1* back=0)
+{
+ TFile* file = TFile::Open(Form("Lego_%s.root", which), "READ");
+ if (!file) {
+ cerr << "Couldn't open the file 'Lego_" << which << ".root"
+ << endl;
+ return 0;
+ }
+ TH1D* h1d = static_cast<TH1D*>(file->Get(Form("h%s_py", what)));
+ if (!h1d) {
+ cerr << "Couldn't find h" << what << "_py in "
+ << "Lego_" << which << ".root" << endl;
+ return 0;
+ }
+
+ TAxis* xaxis = h1d->GetXaxis();
+ Int_t n = xaxis->GetNbins();
+ TArrayF bins(n-1);
+ for (Int_t i = n-1; i > 1; i--) {
+ Float_t ang = xaxis->GetBinUpEdge(i);
+ Float_t eta = deg2eta(ang);
+ bins[n-i-1] = deg2eta(xaxis->GetBinUpEdge(i));
+ }
+ bins[n-2] = deg2eta(xaxis->GetBinLowEdge(2));
+
+ TH1F* heta = new TH1F(Form("%s_eta", what), h1d->GetTitle(),
+ n-2, bins.fArray);
+ heta->SetXTitle("#eta");
+ heta->GetXaxis()->SetTitleSize(0);
+ heta->GetXaxis()->SetTitleOffset(1.5);
+ heta->GetXaxis()->SetTitleFont(132);
+ heta->GetXaxis()->SetLabelFont(132);
+ heta->SetYTitle(Form("%s per degree", h1d->GetTitle()));
+ heta->GetYaxis()->SetTitleOffset(1.5);
+ heta->GetYaxis()->SetTitleFont(132);
+ heta->GetYaxis()->SetLabelFont(132);
+ heta->GetYaxis()->SetTitleSize(0);
+ heta->GetZaxis()->SetTitle(heta->GetTitle());
+ heta->GetZaxis()->SetTitleOffset(1.5);
+ heta->GetZaxis()->SetTitleFont(132);
+ heta->GetZaxis()->SetLabelFont(132);
+ heta->GetYaxis()->SetTitleOffset(1.5);
+ for (Int_t i = 2; i < n; i++) {
+ Float_t ang = xaxis->GetBinUpEdge(i);
+ Float_t eta = deg2eta(ang);
+ Float_t y = h1d->GetBinContent(i);
+ Float_t j = heta->FindBin(eta);
+ if (back) {
+ Float_t b = back->GetBinContent(j);
+ if (y - b <= 0) y = .000001;
+ else y = y - b;
+ }
+ // cout << i << ": " << ang << " -> " << eta << " = " << y << endl;
+ heta->SetBinContent(j, y);
+ }
+ return heta;
+}
+
+void
+drawLego(const char* what="abso")
+{
+ gStyle->SetOptTitle(0);
+ gStyle->SetOptStat(0);
+ gStyle->SetPadLeftMargin(0.15);
+ gStyle->SetPadRightMargin(0.05);
+ gStyle->SetPadTopMargin(0.05);
+ gStyle->SetPadBottomMargin(0.15);
+ gStyle->SetLabelFont(132, "xyz");
+ gStyle->SetTitleFont(132, "xyz");
+ gStyle->SetTitleOffset(1.5, "y");
+
+ TH1* nothing = getHisto("Nothing", what);
+ TH1* its = getHisto("ITS", what, nothing);
+ TH1* fmd = getHisto("FMD", what, nothing);
+ TH1* pipe = getHisto("PIPE", what, nothing);
+ TH1* inner = getHisto("Inner", what);
+
+ if (!inner || !pipe || !fmd || !its || !nothing) {
+ cerr << "Failed to get a histogram!" << endl;
+ return;
+ }
+ TCanvas* c = new TCanvas(Form("single_%s", what),
+ Form("Single %s", what), 800, 800);
+ c->SetFillColor(0);
+ c->SetBorderMode(0);
+ c->SetBorderSize(0);
+ TPad* p1 = new TPad("p1", "p1", 0.0, 0.5, 0.5, 1.0, 0, 0, 0);
+ TPad* p2 = new TPad("p2", "p2", 0.5, 0.5, 1.0, 1.0, 0, 0, 0);
+ TPad* p3 = new TPad("p3", "p3", 0.0, 0.0, 0.5, 0.5, 0, 0, 0);
+ TPad* p4 = new TPad("p4", "p4", 0.5, 0.0, 1.0, 0.5, 0, 0, 0);
+
+ TLatex* latex = new TLatex(0,0,"");
+ latex->SetTextFont(132);
+
+ float logmax = inner->GetMaximum();
+ float logmin = .001 * logmax;
+ float latexy = 6 * logmin;
+
+ c->cd();
+ p1->SetLogy();
+ p1->SetGridy();
+ p1->SetTopMargin(0.15);
+ p1->SetBottomMargin(0);
+ p1->SetRightMargin(0);
+ p1->Draw();
+ p1->cd();
+ fmd->SetFillColor(4);
+ fmd->GetYaxis()->SetRangeUser(logmin, logmax);
+ fmd->GetYaxis()->SetTitleSize(.04);
+ fmd->Draw();
+ latex->DrawLatex(-1, latexy, "FMD only");
+
+ c->cd();
+ p2->SetLogy();
+ p2->SetGridy();
+ p2->SetTopMargin(0.15);
+ p2->SetBottomMargin(0);
+ p2->SetLeftMargin(0);
+ p2->SetRightMargin(0.15);
+ p2->Draw();
+ p2->cd();
+ its->SetFillColor(2);
+ its->GetYaxis()->SetRangeUser(logmin, logmax);
+ its->Draw();
+ latex->DrawLatex(-1, latexy, "ITS only");
+
+ c->cd();
+ p3->SetLogy();
+ p3->SetGridy();
+ p3->SetTopMargin(0);
+ p3->SetLeftMargin(0.15);
+ p3->SetRightMargin(0);
+ p3->Draw();
+ p3->cd();
+ pipe->SetFillColor(3);
+ pipe->GetYaxis()->SetRangeUser(logmin, logmax);
+ pipe->Draw();
+ latex->DrawLatex(-1, latexy, "PIPE only");
+
+ c->cd();
+ p4->SetLogy();
+ p4->SetGridy();
+ p4->SetTopMargin(0);
+ p4->SetLeftMargin(0.);
+ p4->SetRightMargin(0.15);
+ p4->Draw();
+ p4->cd();
+ inner->GetYaxis()->SetRangeUser(logmin, logmax);
+ inner->GetXaxis()->SetTitleSize(.04);
+ inner->SetFillColor(5);
+ inner->Draw();
+ latex->DrawLatex(-1, latexy, "PIPE, ITS, FMD and Air");
+
+ c->Modified();
+ c->cd();
+ c->Print(Form("%s_single.png", what));
+
+ TCanvas* accum = new TCanvas(Form("accum_%s", what),
+ Form("Accumalted %s",what),
+ 800, 500);
+ accum->SetLogy();
+ accum->SetFillColor(0);
+ accum->SetBorderMode(0);
+ accum->SetBorderSize(0);
+
+ THStack* stack = new THStack("stack", "Stack");
+ nothing->SetFillColor(6);
+ stack->Add(nothing);
+ stack->Add(pipe);
+ stack->Add(fmd);
+ stack->Add(its);
+
+ TLegend* legend = new TLegend(.15, .65, .27, .95);
+ legend->SetFillColor(0);
+ legend->SetBorderSize(1);
+ legend->AddEntry(its, "ITS", "f");
+ legend->AddEntry(fmd, "FMD", "f");
+ legend->AddEntry(pipe, "PIPE", "f");
+ legend->AddEntry(nothing, "Air", "f");
+
+ stack->SetMinimum(nothing->GetMinimum());
+ stack->SetMaximum(logmax);
+ stack->Draw();
+ stack->GetXaxis()->SetTitle("#eta");
+ // stack->GetYaxis()->SetRangeUser(, );
+ stack->GetYaxis()->SetTitle(fmd->GetTitle());
+ legend->Draw();
+ accum->Modified();
+
+ accum->Modified();
+ accum->cd();
+ // accum->Print(Form("%s_accum.eps", what));
+ accum->Print(Form("%s_accum.png", what));
+
+}
+
+
+
+void
+DrawLego()
+{
+ drawLego("abso");
+ drawLego("radl");
+ drawLego("gcm2");
+}
+
+
--- /dev/null
+//
+// Script to make ROOT files with LEGO information histograms
+//
+TH1*
+process(TFile* file, const char* name, const char* opt)
+{
+ TH2F* h2d = static_cast<TH2F*>(file->Get(name));
+ if (!h2d) {
+ Error("process", "Couldn't get %s from %s", name, file->GetName());
+ return 0;
+ }
+ TH1D* h1d = h2d->ProjectionY();
+ h1d->SetTitle(Form("%s", h1d->GetTitle()));
+ h1d->Scale(1. / 360.);
+
+ return h1d;
+
+ return heta;
+}
+
+void
+MakeLego(const Char_t* what)
+{
+ TString config("FMD/scripts/ConfigInner.C");
+ TString opt(what);
+ if (opt == "ITS")
+ config = "FMD/scripts/ConfigItsOnly.C";
+ else if (opt == "PIPE")
+ config = "FMD/scripts/ConfigPipeOnly.C";
+ else if (opt == "FMD")
+ config = "FMD/scripts/ConfigFmdOnly.C";
+ else if (opt == "Nothing")
+ config = "FMD/scripts/ConfigNothing.C";
+ else
+ opt = "Inner";
+
+
+ cout << "Running AliRun::RunLego(" << config
+ << ",180,0,180,360,0,360,0,10000,0,10000); " << endl;
+
+ gAlice->RunLego(config.Data(), 180, 0, 180, 360, 0, 360, 0, 100000,
+ 100000000, 0);
+
+ TFile* galice = TFile::Open("galice.root", "READ");
+ TFile* output = TFile::Open(Form("Lego_%s.root", opt.Data()),"RECREATE");
+ TH1F* habso_eta = process(galice, "habso", opt.Data());
+ TH1F* hradl_eta = process(galice, "hradl", opt.Data());
+ TH1F* hgcm2_eta = process(galice, "hgcm2", opt.Data());
+ hgcm2_eta->Draw();
+
+ output->Write();
+ output->Close();
+ galice->Close();
+}
--- /dev/null
+//
+// $Id$
+//
+// Script to make a class derived from AliFMDMap.
+//
+#ifndef __CINT__
+#include <fstream>
+#include <TString.h>
+#include <TDatime.h>
+#include <TSystem.h>
+#include <iostream>
+using namespace std;
+#endif
+
+void MakeMap(const Char_t* type="Int_t", const Char_t* name=0)
+{
+ TString base;
+ TString ttype(type);
+ if (ttype.EndsWith("_t")) {
+ Ssiz_t undert = ttype.Index("_t");
+ ttype.Remove(undert);
+ }
+ if (!name)
+ base = Form("AliFMD%sMap", ttype.Data());
+ else
+ base = name;
+
+ cout << "Base name is " << base << endl;
+
+ TString decl_name(Form("%s.h", base.Data()));
+ TString impl_name(Form("%s.cxx", base.Data()));
+ ofstream decl(decl_name.Data());
+ ofstream impl(impl_name.Data());
+
+ if (!decl) {
+ cerr << "Cannot open declaration file " << decl_name << endl;
+ return;
+ }
+ if (!impl) {
+ cerr << "Cannot open implementation file " << impl_name << endl;
+ return;
+ }
+
+ TDatime now;
+ cout << "The time is now " << now.AsString() << endl;
+ UserGroup_t* uinfo = gSystem->GetUserInfo();
+ TString uname(uinfo->fRealName);
+ Ssiz_t comma = uname.Index(",");
+ if (comma != kNPOS) uname.Remove(comma);
+ cout << "User's name is " << uname << endl;
+ TString guard(base);
+ guard.Append("_h");
+ guard.ToUpper();
+
+ cout << "Writing declaration file " << decl_name << " ... "
+ << flush;
+ decl << "#ifndef " << guard << "\n";
+ decl << "#define " << guard << "\n";
+ decl << "/* Copyright (c) " << now.GetYear() ;
+ decl << ", ALICE Experiment @ CERN.\n" ;
+ decl << " * All rights reserved\n";
+ decl << " * See " << impl_name << " for full copyright notice\n";
+ decl << " * \n" ;
+ decl << " * Created " << now.AsString() << " by " << uname << "\n";
+ decl << " */\n";
+ decl << "/* $Id$ */\n";
+ decl << "//__________________________________________________________\n";
+ decl << "// \n";
+ decl << "// Map of " << type << " for each FMD strip\n" ;
+ decl << "// \n";
+ decl << "#ifndef ALIFMDMAP_H\n";
+ decl << "# include <AliFMDMap.h>\n";
+ decl << "#endif\n\n";
+ decl << "class " << base << " : public AliFMDMap\n";
+ decl << "{\n";
+ decl << "public:\n";
+ decl << " " << base << "(const " << base << "& other);\n";
+ decl << " " << base << "(size_t maxDet = kMaxDetectors,\n";
+ decl << " size_t maxRing = kMaxRings,\n";
+ decl << " size_t maxSec = kMaxSectors,\n";
+ decl << " size_t maxStr = kMaxStrips);\n";
+ decl << " virtual ~" << base << "() { delete [] fData; }\n";
+ decl << " " << base << "& operator=(const " << base << "& other);\n";
+ decl << " virtual void Clear(const " << type << "& v=" << type << "());\n";
+ decl << " virtual " << type << "& operator()(UShort_t det,\n";
+ decl << " Char_t ring,\n";
+ decl << " UShort_t sec,\n";
+ decl << " UShort_t str);\n";
+ decl << " virtual const " << type << "& operator()(UShort_t det,\n";
+ decl << " Char_t ring,\n";
+ decl << " UShort_t sec,\n";
+ decl << " UShort_t str) const;\n";
+ decl << "protected:\n";
+ decl << " " << type << "* fData; // The Data\n";
+ decl << " ClassDef(" << base << ",1) // Map of " << type ;
+ decl << " data per strip\n" ;
+ decl << "};\n\n";
+ decl << "#endif\n";
+ decl << "//__________________________________________________________\n";
+ decl << "// \n";
+ decl << "// Local Variables:\n";
+ decl << "// mode: C++\n";
+ decl << "// End:\n";
+ decl << "//" << endl;;
+ decl.close();
+ cout << "done" << endl;
+
+ cout << "Writing implementation file " << impl_name << " ... "
+ << flush;
+ impl << "/**************************************************************\n";
+ impl << " * Copyright(c) 1998-1999, ALICE Experiment at CERN. *\n";
+ impl << " * All rights reserved. *\n";
+ impl << " * *\n";
+ impl << " * Author: The ALICE Off-line Project. *\n";
+ impl << " * Contributors are mentioned in the code where appropriate. *\n";
+ impl << " * *\n";
+ impl << " * Permission to use, copy, modify and distribute this *\n";
+ impl << " * software and its documentation strictly for non-commercial *\n";
+ impl << " * purposes is hereby granted without fee, provided that the *\n";
+ impl << " * above copyright notice appears in all copies and that both *\n";
+ impl << " * the copyright notice and this permission notice appear in *\n";
+ impl << " * the supporting documentation. The authors make no claims *\n";
+ impl << " * about the suitability of this software for any purpose. It *\n";
+ impl << " * is provided \"as is\" without express or implied warranty. *\n";
+ impl << " **************************************************************/\n";
+ impl << "/* $Id$ */\n";
+ impl << "//__________________________________________________________\n";
+ impl << "// \n";
+ impl << "// Map of per strip " << type << " information\n";
+ impl << "// \n";
+ impl << "// Created " << now.AsString() << " by " << uname << "\n";
+ impl << "// \n";
+ impl << "#include \"" << decl_name << "\"\t//" << guard << "\n";
+ impl << "//__________________________________________________________\n";
+ impl << "ClassImp(" << base << ");\n";
+ impl << "//__________________________________________________________\n";
+ impl << base << "::" << base << "(const " << base << "& other)\n";
+ impl << " : AliFMDMap(other.fMaxDetectors,\n";
+ impl << " other.fMaxRings,\n";
+ impl << " other.fMaxSectors,\n";
+ impl << " other.fMaxStrips),\n";
+ impl << " fData(0)\n";
+ impl << "{\n";
+ impl << " // Copy constructor\n";
+ impl << " fData = new " << type << "[fMaxDetectors * fMaxRings ";
+ impl << "* fMaxSectors * fMaxStrips];\n" ;
+ impl << " for (size_t i = 0; i < fMaxDetectors * fMaxRings ";
+ impl << "* fMaxSectors * fMaxStrips; i++)\n";
+ impl << " fData[i] = other.fData[i];\n";
+ impl << "}\n\n";
+ impl << "//__________________________________________________________\n";
+ impl << base << "::" << base << "(size_t maxDet,\n";
+ impl << " size_t maxRing,\n";
+ impl << " size_t maxSec,\n";
+ impl << " size_t maxStr)\n";
+ impl << " : AliFMDMap(maxDet, maxRing, maxSec, maxStr),\n";
+ impl << " fData(0)\n";
+ impl << "{\n";
+ impl << " // Constructor.\n";
+ impl << " // Parameters:\n";
+ impl << " //\tmaxDet\tMaximum number of detectors\n";
+ impl << " //\tmaxRing\tMaximum number of rings per detector\n";
+ impl << " //\tmaxSec\tMaximum number of sectors per ring\n";
+ impl << " //\tmaxStr\tMaximum number of strips per sector\n";
+ impl << " fData = new " << type << "[fMaxDetectors * fMaxRings ";
+ impl << "* fMaxSectors * fMaxStrips];\n" ;
+ impl << " Clear();\n";
+ impl << "}\n\n";
+ impl << "//__________________________________________________________\n";
+ impl << base << "&\n";
+ impl << base << "::operator=(const " << base << "& other)\n";
+ impl << "{\n";
+ impl << " // Assignment operator \n";
+ impl << " fMaxDetectors = other.fMaxDetectors;\n";
+ impl << " fMaxRings = other.fMaxRings;\n";
+ impl << " fMaxSectors = other.fMaxSectors;\n";
+ impl << " fMaxStrips = other.fMaxStrips;\n";
+ impl << " if (fData) delete [] fData;\n";
+ impl << " fData = new " << type << "[fMaxDetectors * fMaxRings ";
+ impl << "* fMaxSectors * fMaxStrips];\n" ;
+ impl << " for (size_t i = 0; i < fMaxDetectors * fMaxRings ";
+ impl << "* fMaxSectors * fMaxStrips; i++)\n";
+ impl << " fData[i] = other.fData[i];\n";
+ impl << "}\n\n" ;
+ impl << "//__________________________________________________________\n";
+ impl << "void\n";
+ impl << base << "::Clear(const " << type << "& val)\n";
+ impl << "{\n";
+ impl << " // Reset map to val\n";
+ impl << " for (size_t i = 0; i < fMaxDetectors * fMaxRings ";
+ impl << "* fMaxSectors * fMaxStrips; i++)\n";
+ impl << " fData[i] = val;\n";
+ impl << "}\n\n" ;
+ impl << "//__________________________________________________________\n";
+ impl << type << "&\n";
+ impl << base << "::operator()(UShort_t det, Char_t ring, UShort_t sec, " ;
+ impl << "UShort_t str)\n" ;
+ impl << "{\n" ;
+ impl << " // Get data\n";
+ impl << " // Parameters:\n";
+ impl << " //\tdet\tDetector #\n";
+ impl << " //\tring\tRing ID\n";
+ impl << " //\tsec\tSector #\n";
+ impl << " //\tstr\tStrip #\n" ;
+ impl << " // Returns appropriate data\n";
+ impl << " return fData[CalcIndex(det, ring, sec, str)];\n";
+ impl << "}\n\n";
+ impl << "//__________________________________________________________\n";
+ impl << "const " << type << "&\n";
+ impl << base << "::operator()(UShort_t det, Char_t ring, UShort_t sec, " ;
+ impl << "UShort_t str) const\n" ;
+ impl << "{\n" ;
+ impl << " // Get data\n";
+ impl << " // Parameters:\n";
+ impl << " //\tdet\tDetector #\n";
+ impl << " //\tring\tRing ID\n";
+ impl << " //\tsec\tSector #\n";
+ impl << " //\tstr\tStrip #\n" ;
+ impl << " // Returns appropriate data\n";
+ impl << " return fData[CalcIndex(det, ring, sec, str)];\n";
+ impl << "}\n\n";
+ impl << "//__________________________________________________________\n";
+ impl << "// \n";
+ impl << "// EOF\n";
+ impl << "// \n";
+ impl << endl;
+ impl.close();
+ cout << "done" << endl;
+}
+
+
+
+
+#ifndef __CINT__
+int main()
+{
+ makemap();
+ return 0;
+}
+#endif
--- /dev/null
+//
+// Small script to test consistency of writing and reading raw data.
+//
+void
+RawTest()
+{
+ gRandom->SetSeed(12345);
+ Int_t sampleRate = 3;
+ Int_t channelWidth = 128;
+ Float_t shapingTime = 5;
+ UInt_t maxAdc = (1 << 10);
+ UInt_t threshold = (1 << 8);
+ TArrayI outData(sampleRate * channelWidth);
+
+ Float_t lastTotalCharge = 0;
+ Int_t ok = 0;
+ for (Int_t channel = 0; channel < channelWidth; channel++) {
+ Float_t totalCharge = gRandom->Uniform(0, 1);
+
+ for (Int_t sample = 0; sample < sampleRate; sample++) {
+ Float_t time = Float_t(sample) / sampleRate;
+ Float_t charge = (totalCharge + (lastTotalCharge - totalCharge)
+ * TMath::Exp(-shapingTime * time));
+ UInt_t adc = UInt_t(maxAdc * charge);
+ outData[channel * sampleRate + sample] = adc;
+ if (adc > threshold) ok++;
+ }
+ lastTotalCharge = totalCharge;
+ }
+ std::cout << "Total of " << outData.fN << " samples of which "
+ << ok << " of them are above threshold (" << threshold
+ << ")" << std::endl;
+
+ {
+ AliAltroBuffer buffer("FMD_4096.ddl", 1);
+ buffer.WriteDataHeader(kTRUE, kFALSE);
+ buffer.WriteChannel(0, 0, 0, outData.fN, outData.fArray, threshold);
+ buffer.Flush();
+ buffer.WriteDataHeader(kFALSE, kFALSE);
+ }
+
+ AliRawReader* reader = new AliRawReaderFile(-1);
+ if (!reader) {
+ std::cerr << "Failed to make AliRawReader" << endl;
+ return 0;
+ }
+ AliFMDRawStream input(reader, sampleRate);
+ reader->Select(AliFMD::kBaseDDL >> 8);
+
+ Int_t oldDDL = -1;
+ Int_t count = 0;
+ UShort_t detector = 1; // Must be one here
+ UShort_t oldDetector = 0;
+ // Loop over data in file
+ Bool_t next = kTRUE;
+
+ // local Cache
+ TArrayI counts(10);
+ counts.Reset(-1);
+ Int_t offset = 0;
+
+ TArrayI inputData(sampleRate * channelWidth);
+ while (next) {
+ next = input.Next();
+
+ if (!next) break;
+
+ Int_t channel = input.Strip();
+ Int_t sample = input.Sample();
+ inputData[channel * sampleRate + sample] = input.Count();
+ count++;
+
+ Int_t in = inputData[channel * sampleRate + sample];
+ Int_t out = outData[channel * sampleRate + sample];
+ std::cout << "[\t" << channel << ",\t" << sample << "]\t"
+ << out << "\t" << in << std::flush;
+ if (out >= threshold && in != out) std::cout << "\tBad" << std::flush;
+ std::cout << std::endl;
+ }
+
+ std::cout << "Read " << count << " values" << std::endl;
+#if 1
+ for (Int_t channel = channelWidth - 1; channel > 0; channel--) {
+ for (Int_t sample = sampleRate - 1; sample > 0; sample--) {
+ Int_t in = inputData[channel * sampleRate + sample];
+ Int_t out = outData[channel * sampleRate + sample];
+ std::cout << "[\t" << channel << ",\t" << sample << "]\t"
+ << out << "\t" << in << std::flush;
+ if (out >= threshold && in != out) std::cout << "\tBad" << std::flush;
+ std::cout << std::endl;
+ }
+ }
+#endif
+}
+
+
+
+
+
+
+
+
+
--- /dev/null
+//
+// Script to digit multiplicity information to std::cout.
+//
+void
+ShowDigits()
+{
+ AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
+ runLoader->LoadgAlice();
+ runLoader->LoadHeader();
+ gAlice = runLoader->GetAliRun();
+ AliFMD* fmd = static_cast<AliFMD*>(gAlice->GetDetector("FMD"));
+ AliLoader* fmdLoader = runLoader->GetLoader("FMDLoader");
+ fmdLoader->LoadDigits("READ");
+
+ Int_t nEvents = runLoader->TreeE()->GetEntries();
+ for (Int_t event = 0; event < nEvents; event++) {
+ cout << "Event # " << event << endl;
+ runLoader->GetEvent(event);
+ TClonesArray* digits = 0;
+ TTree* treeD = fmdLoader->TreeD();
+ TBranch* branch = treeD->GetBranch("FMD");
+ branch->SetAddress(&digits);
+
+ Int_t total = 0;
+ Int_t nEntries = treeD->GetEntries();
+ for (Int_t entry = 0; entry < nEntries; entry++) {
+ cout << " Entry # " << entry << endl;
+ treeD->GetEntry(entry);
+
+ Int_t nDigits = digits->GetLast();
+ for (Int_t i = 0; i < nDigits; i++) {
+ // cout << " Digit # " << i << endl;
+ AliFMDDigit* digit = static_cast<AliFMDDigit*>(digits->UncheckedAt(i));
+ if (digit->Counts() > 12) {
+ digit->Print();
+ total++;
+ }
+ }
+ }
+ cout << "Total number of digits: " << total << endl;
+ }
+}
--- /dev/null
+//
+// Script to dump hit information to std::cout.
+//
+ShowHits()
+{
+ AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
+ runLoader->LoadgAlice();
+ runLoader->LoadHeader();
+ gAlice = runLoader->GetAliRun();
+ AliFMD* fmd = static_cast<AliFMD*>(gAlice->GetDetector("FMD"));
+ AliLoader* fmdLoader = runLoader->GetLoader("FMDLoader");
+ fmdLoader->LoadHits("READ");
+
+ Int_t nEvents = runLoader->TreeE()->GetEntries();
+ for (Int_t event = 0; event < nEvents; event++) {
+ cout << "Event # " << event << endl;
+ runLoader->GetEvent(event);
+ TClonesArray* hits = 0;
+ TTree* treeH = fmdLoader->TreeH();
+ TBranch* branch = treeH->GetBranch("FMD");
+ branch->SetAddress(&hits);
+
+ Int_t total = 0;
+ Int_t nEntries = treeH->GetEntries();
+ for (Int_t entry = 0; entry < nEntries; entry++) {
+ // cout << " Entry # " << entry << endl;
+ treeH->GetEntry(entry);
+
+ Int_t nHits = hits->GetEntries();
+ for (Int_t i = 0; i < nHits; i++) {
+ cout << " Hit # " << i << "/" << nHits << "\t" << flush;
+ AliFMDHit* hit = static_cast<AliFMDHit*>(hits->UncheckedAt(i));
+ hit->Print();
+ total++;
+ }
+ }
+ cout << "Total number of hits: " << total << endl;
+ }
+}
--- /dev/null
+//
+// Script to dump multiplicity information to std::cout.
+//
+void
+ShowMult()
+{
+ AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
+ runLoader->LoadgAlice();
+ runLoader->LoadHeader();
+ gAlice = runLoader->GetAliRun();
+ AliFMD* fmd = static_cast<AliFMD*>(gAlice->GetDetector("FMD"));
+ AliLoader* fmdLoader = runLoader->GetLoader("FMDLoader");
+ fmdLoader->LoadRecPoints("READ");
+
+ Int_t nEvents = runLoader->TreeE()->GetEntries();
+ for (Int_t event = 0; event < nEvents; event++) {
+ cout << "Event # " << event << endl;
+ runLoader->GetEvent(event);
+ TClonesArray* multStrips = 0;
+ TClonesArray* multRegions = 0;
+ TTree* treeR = fmdLoader->TreeR();
+ TBranch* branchRegions = treeR->GetBranch("FMDPoisson");
+ TBranch* branchStrips = treeR->GetBranch("FMDNaiive");
+ branchRegions->SetAddress(&multRegions);
+ branchStrips->SetAddress(&multStrips);
+
+ Int_t total = 0;
+ Int_t nEntries = treeR->GetEntries();
+ for (Int_t entry = 0; entry < nEntries; entry++) {
+ cout << " Entry # " << entry << endl;
+ treeR->GetEntry(entry);
+
+ Int_t nMults = multStrips->GetLast();
+ for (Int_t i = 0; i < nMults; i++) {
+ // cout << " Digit # " << i << endl;
+ AliFMDMultStrip* mult =
+ static_cast<AliFMDMultStrip*>(multStrips->UncheckedAt(i));
+ if (mult->Particles() > 0) mult->Print();
+ }
+
+ nMults = multRegions->GetLast();
+ for (Int_t i = 0; i < nMults; i++) {
+ // cout << " Digit # " << i << endl;
+ AliFMDMultStrip* mult =
+ static_cast<AliFMDMultStrip*>(multRegions->UncheckedAt(i));
+ if (mult->Particles() > 0) mult->Print();
+ }
+ }
+ }
+}
--- /dev/null
+//
+// Script to read a raw data file, and dump it to std::cout
+//
+#include <iomanip>
+
+void
+ShowRaw(Int_t det=2, bool verbose=false, Int_t event=0)
+{
+ TString file(Form("raw%d/FMD_%d.ddl", event, AliFMD::kBaseDDL + det - 1));
+
+ std::cout << "Reading raw data file " << file << std::endl;
+
+ TH1* h = new TH1F("rawData", "Raw Data", 90, 0, 90);
+
+
+ // This method creates a text file containing the same information
+ // stored in an Altro file. The information in the text file is
+ // organized pad by pad and and for each pad it consists in a
+ // sequence of bunches (Bunch length +2, Time bin of the last
+ // amplitude sample in the bunch, amplitude values) It is used
+ // mainly //for debugging
+
+ AliAltroBuffer buff(file.Data(),0);
+ Int_t numWords,padNum,rowNum,secNum=0;
+ Int_t value = 0;
+ Int_t zero = 0;
+ // if (!buff.ReadDataHeader()) {
+ // std::cout<< file << " isn't a valid data file!" << std::endl;
+ // }
+
+ while(buff.ReadTrailerBackward(numWords,padNum,rowNum,secNum)){
+ if (verbose)
+ std::cout << "Ring: " << (secNum == 0 ? 'I' : 'O')
+ << " Sector: " << std::setw(2) << rowNum
+ << " Strip: " << std::setw(3) << padNum
+ << " Words: " << std::setw(4) << numWords << std::endl;
+ if (numWords == 0) zero++;
+ if (numWords % 4){
+ if (verbose)
+ std::cout << "Skipping trailer of "
+ << (4 - numWords % 4) << " words" << std::endl;
+ for(Int_t j = 0; j < (4 - numWords % 4); j++)
+ value=buff.GetNextBackWord();
+ }//end if
+ for(Int_t i = 0; i <numWords; i++) {
+ value=buff.GetNextBackWord();
+ if (verbose) {
+ std::cout << std::setw(5) << value << std::flush;
+ if (i % 16 == 15) std::cout << std::endl;
+ }
+ h->Fill(value);
+ }//end for
+ if (verbose)
+ std::cout << std::endl;
+ if (zero > 1) {
+ std::cout << "Error: Read zero channels - should not happen"
+ << std::endl;
+ break;
+ }
+ }//end while
+ h->Draw();
+ return;
+}
--- /dev/null
+//____________________________________________________________________
+//
+// Script I used for rapid prototyping of the FMD3 geometry - in
+// particular the support cone
+//
+//____________________________________________________________________
+TObjArray*
+waferParameters(double dl, double dh, double theta, double r,
+ Bool_t verbose=kFALSE)
+{
+ double tan_theta = tan(theta * TMath::Pi() / 180.);
+ double tan_theta2 = pow(tan_theta,2);
+ double r2 = pow(r,2);
+ double y_A = tan_theta * dl;
+ double x_D = dl + sqrt(r2 - tan_theta2 * pow(dl,2));
+ double x_D_2 = dl - sqrt(r2 - tan_theta2 * pow(dl,2));
+
+ double y_B = sqrt(r2 - pow(dh,2) + 2 * dh * x_D - pow(x_D,2));
+ double x_C = (x_D + sqrt(-tan_theta2 * pow(x_D,2) + r2
+ + r2 * tan_theta2)) / (1 + tan_theta2);
+ double y_C = tan_theta * x_C;
+
+ if (verbose) {
+ cout << "A: (" << dl << "," << y_A << ")" << endl;
+ cout << "B: (" << dh << "," << y_B << ")" << endl;
+ cout << "C: (" << x_C << "," << y_C << ")" << endl;
+ cout << "D: (" << x_D << ",0)" << endl;
+
+ cout << "Recentred at D:" << endl;
+ cout << "A: (" << dl - x_D << "," << y_A << ")" << endl;
+ cout << "B: (" << dh - x_D << "," << y_B << ")" << endl;
+ cout << "C: (" << x_C - x_D << "," << y_C << ")" << endl;
+ }
+
+ TObjArray* verticies = new TObjArray(6);
+ verticies->AddAt(new TVector2(dl, y_A), 0);
+ verticies->AddAt(new TVector2(x_C, y_C), 1);
+ verticies->AddAt(new TVector2(dh, y_B), 2);
+ verticies->AddAt(new TVector2(dh, -y_B), 3);
+ verticies->AddAt(new TVector2(x_C, -y_C), 4);
+ verticies->AddAt(new TVector2(dl, -y_A), 5);
+
+ return verticies;
+}
+
+//____________________________________________________________________
+TShape*
+createModuleShape(const Char_t* name, double rl, double rh, double th,
+ double r, double dz)
+{
+ std::cout << "Creating Module shape for " << name << std::flush;
+ // TShape* virtualShape = new TTUBS(Form("%sVirtual", name),
+ // Form("Virtual %s", name),
+ // "", rl, rh, 1, -th, th);
+
+ TObjArray* v = waferParameters(rl, rh, th, r);
+ TXTRU* moduleShape = new TXTRU(Form("%sModule", name),
+ Form("Module %s", name),
+ "", 6, 2);
+ for (Int_t i = 0; i < 6; i++) {
+ std::cout << "." << std::flush;
+ TVector2* vv = static_cast<TVector2*>(v->At(i));
+ moduleShape->DefineVertex(i, vv->X(), vv->Y());
+ }
+ moduleShape->DefineSection(0, -dz, 1, 0, 0);
+ moduleShape->DefineSection(1, dz, 1, 0, 0);
+ std::cout << std::endl;
+
+ return (TShape*)moduleShape;
+}
+
+//____________________________________________________________________
+TNode*
+createRing(const char* name, double rl, double rh, double th,
+ double siThick, double waferR, double staggering, double z)
+{
+ std::cout << "Creating Ring node for " << name << std::flush;
+ TShape* ringShape = new TTUBE(Form("%sShape", name), "Ring Shape",
+ "", rl, rh, staggering + siThick);
+ TNode* ringNode = new TNode(Form("%sNode", name), "Ring Node",
+ ringShape, 0, 0, z, 0);
+ TShape* moduleShape = createModuleShape(name, rl, rh, th, waferR, siThick);
+ Int_t n = 360 / 2 / th;
+ for (Int_t i = 0; i < n; i++) {
+ std::cout << "." << std::flush;
+ ringNode->cd();
+ Double_t theta = 2 * th * i;
+ Double_t z = (i % 2 ? 0 : staggering);
+ TRotMatrix* rot = new TRotMatrix(Form("%sRotation%02d", name, i),
+ "Rotation", 90, theta, 90,
+ fmod(90 + theta, 360), 0, 0);
+ TNode* moduleNode = new TNode(Form("%sModule%02d", name, i),
+ "Module", moduleShape, 0, 0, z,
+ rot);
+ moduleNode->SetFillColor(5);
+ moduleNode->SetLineColor(5);
+ moduleNode->SetLineWidth(2);
+ }
+ std::cout << std::endl;
+ ringNode->SetVisibility(0);
+ return ringNode;
+}
+
+//____________________________________________________________________
+TNode*
+createSupport(double noseRl, double noseRh, double noseDz, double noseZ,
+ double backRl, double backRh, double backDz, double coneL,
+ double beamW, double beamDz, double flangeR)
+{
+ TShape* noseShape = new TTUBE("noseShape", "Nose Shape", "",
+ noseRl, noseRh, noseDz);
+ TNode* noseNode = new TNode("noseNode", "noseNode", noseShape,
+ 0, 0, noseZ - noseDz, 0);
+ noseNode->SetLineColor(0);
+
+ Double_t zdist = coneL - 2 * backDz - 2 * noseDz;
+ Double_t tdist = backRh - noseRh;
+ Double_t beamL = TMath::Sqrt(zdist * zdist + tdist * tdist);
+ Double_t theta = -TMath::ATan2(tdist, zdist);
+
+
+ TShape* backShape = new TTUBE("backShape", "Back Shape", "",
+ backRl, backRh, backDz);
+ TNode* backNode = new TNode("backNode", "backNode", backShape,
+ 0, 0, noseZ - 2 * noseDz - zdist - backDz, 0);
+ backNode->SetLineColor(0);
+
+
+ TShape* beamShape = new TBRIK("beamShape", "beamShape", "",
+ beamDz, beamW / 2 , beamL / 2);
+ Int_t n = 8;
+ Double_t r = noseRl + tdist / 2;
+ for (Int_t i = 0; i < n; i++) {
+ Double_t phi = 360. / n * i;
+ Double_t t = 180. * theta / TMath::Pi();
+ TRotMatrix* beamRotation = new TRotMatrix(Form("beamRotation%d", i),
+ Form("beamRotation%d", i),
+ 180-t, phi, 90, 90+phi,
+ t, phi);
+ TNode* beamNode = new TNode(Form("beamNode%d", i),
+ Form("beamNode%d", i), beamShape,
+ r * TMath::Cos(phi / 180 * TMath::Pi()),
+ r * TMath::Sin(phi / 180 * TMath::Pi()),
+ noseZ - 2 * noseDz - zdist / 2, beamRotation);
+ beamNode->SetLineColor(0);
+ }
+
+ Double_t flangel = (flangeR - backRh) / 2;
+ TShape* flangeShape = new TBRIK("flangeShape", "FlangeShape", "",
+ flangel, beamW / 2, backDz);
+ n = 4;
+ r = backRh + flangel;
+ for (Int_t i = 0; i < n; i++) {
+ Double_t phi = 360. / n * i + 180. / n;
+ TRotMatrix* flangeRotation = new TRotMatrix(Form("flangeRotation%d", i),
+ Form("Flange Rotation %d", i),
+ 90, phi, 90, 90+phi, 0, 0);
+ TNode* flangeNode = new TNode(Form("flangeNode%d", i),
+ Form("flangeNode%d", i),
+ flangeShape,
+ r * TMath::Cos(phi / 180 * TMath::Pi()),
+ r * TMath::Sin(phi / 180 * TMath::Pi()),
+ noseZ - 2 * noseDz - zdist - backDz,
+ flangeRotation);
+ flangeNode->SetLineColor(0);
+
+ }
+ return 0;
+}
+
+
+
+
+//____________________________________________________________________
+void
+SimpleGeometry()
+{
+ TGeometry* geometry = new TGeometry("geometry","geometry");
+ TShape* topShape = new TBRIK("topShape", "topShape", "", 40, 40, 150);
+ TNode* topNode = new TNode("topNode", "topNode", topShape, 0, 0, 0, 0);
+ topNode->SetVisibility(0);
+ topNode->cd();
+
+ Double_t waferR = 13.4 / 2;
+ Double_t siThick = .03;
+ Double_t staggering = 1;
+ Double_t innerRh = 17.2;
+ Double_t innerRl = 4.3;
+ Double_t innerTh = 18;
+ Double_t innerZ = -62.8;
+ Double_t outerRh = 28;
+ Double_t outerRl = 15.6;
+ Double_t outerTh = 9;
+ Double_t outerZ = -75.2;
+ Double_t noseRl = 5.5;
+ Double_t noseRh = 6.7;
+ Double_t noseDz = 2.8 / 2;
+ Double_t noseZ = -46;
+ Double_t coneL = 30.9;
+ Double_t backRl = 61 / 2;
+ Double_t backRh = 66.8 /2;
+ Double_t backDz = 1.4 / 2;
+ Double_t beamDz = .5 / 2;
+ Double_t beamW = 6;
+ Double_t flangeR = 49.25;
+
+ Double_t zdist = coneL - 2 * backDz - 2 * noseDz;
+ Double_t tdist = backRh - noseRh;
+ Double_t alpha = tdist / zdist;
+
+ Double_t x, rl, rh, z;
+ z = noseZ - coneL / 2;
+ TPCON* fmd3Shape = new TPCON("fmd3Shape", "FMD 3 Shape", "", 0, 360, 7);
+
+ x = noseZ;
+ rl = noseRl;
+ rh = noseRh;
+ fmd3Shape->DefineSection(0, x - z, rl, rh);
+
+ x = noseZ-2*noseDz;
+ fmd3Shape->DefineSection(1, x - z, rl, rh);
+
+ x = innerZ - staggering - siThick;
+ rl = innerRl;
+ rh = noseRh + alpha * TMath::Abs(x-noseZ + 2 * noseDz);
+ fmd3Shape->DefineSection(2, x - z, rl, rh);
+
+ x = outerZ;
+ rl = outerRl;
+ rh = backRh;
+ fmd3Shape->DefineSection(3, x - z, rl, rh);
+
+ x = noseZ - zdist - 2 * noseDz;
+ rl = outerRl;
+ rh = backRh;
+ fmd3Shape->DefineSection(4, x - z, rl, rh);
+
+ x = noseZ - zdist - 2 * noseDz;
+ rl = outerRl;
+ rh = flangeR;
+ fmd3Shape->DefineSection(5, x - z, rl, rh);
+
+ x = noseZ - coneL;
+ rl = outerRl;
+ rh = flangeR;
+ fmd3Shape->DefineSection(6, x - z, rl, rh);
+
+ TNode* fmd3Node = new TNode("fmd3Node", "FMD3 Node", fmd3Shape,
+ 0, 0, z, 0);
+ fmd3Node->SetLineColor(3);
+ fmd3Node->SetVisibility(1);
+
+ fmd3Node->cd();
+ TNode* innerNode = createRing("inner", innerRl, innerRh, innerTh,
+ siThick, waferR, staggering, innerZ-z);
+
+
+ fmd3Node->cd();
+ TNode* outerNode = createRing("outer", outerRl, outerRh, outerTh,
+ siThick, waferR, staggering, outerZ-z);
+
+
+ fmd3Node->cd();
+ TNode* supportNode = createSupport(noseRl, noseRh, noseDz, noseZ-z,
+ backRl, backRh, backDz, coneL,
+ beamW, beamDz, flangeR);
+
+ TCanvas* c = new TCanvas("c", "c", 800, 800);
+ c->SetFillColor(1);
+ geometry->Draw();
+ // c->x3d("ogl");
+}
--- /dev/null
+//
+// Script to try to fit the reponse function of the VA1 signals, based
+// on a finite number of ALTRO samples.
+//
+void
+VA1Response(Int_t n=2, Float_t B=5, Float_t dc=.01, Bool_t errors=kFALSE)
+{
+
+ TF1* response = new TF1("response", "[0] * (1 - exp(-[1] * x))", 0, 1.4);
+ response->SetParameters(1, B);
+ response->SetParNames("A", "B");
+ response->SetLineColor(2);
+
+ TF1* fit = new TF1("fit", "[0] * (1 - exp(-[1] * x))", 0, 1);
+ fit->SetParameters(.5, B/2);
+ fit->SetParNames("A", "B");
+ fit->SetLineColor(3);
+
+ TGraph* graph = 0;
+ if (errors) graph = new TGraphErrors(n);
+ else graph = new TGraph(n);
+ for (Int_t i = 0; i < n; i++) {
+ Float_t t = Float_t(i + 1) / n;
+ Float_t c = gRandom->Gaus(response->Eval(t), dc);
+ graph->SetPoint(i, t, c);
+ if (errors) ((TGraphErrors*)graph)->SetPointError(i, 0, dc);
+ }
+
+ response->Draw();
+ response->GetHistogram()->GetYaxis()->SetRangeUser(0, 1.4);
+ response->GetHistogram()->GetXaxis()->SetRangeUser(0, 1.4);
+ graph->Draw("P*");
+ TString fitOpt("E");
+ if (!errors) fitOpt.Append("W");
+ graph->Fit("fit", fitOpt.Data());
+ graph->Fit("fit", fitOpt.Data());
+
+ std::cout << "Chi^2/NDF = " << fit->GetChisquare() << "/" << fit->GetNDF()
+ << " = " << std::flush;
+ if (fit->GetNDF() == 0)
+ std::cout << " undefined!" << std::endl;
+ else
+ std::cout << (fit->GetChisquare() / fit->GetNDF()) << std::endl;
+ std::cout << "f(t) = "
+ << fit->GetParameter(0) << "+/-" << fit->GetParError(0)
+ << " * (1 - exp("
+ << fit->GetParameter(1) << "+/-" << fit->GetParError(1)
+ << " * t))" << std::endl;
+}
--- /dev/null
+//
+// Small script that shows a signal train from a VA1 pre-amp.
+//
+void
+VA1Train()
+{
+ gStyle->SetOptTitle(0);
+ gStyle->SetOptFit(0);
+ gStyle->SetOptStat(0);
+
+ TCanvas* c = new TCanvas("c", "C", 800, 400);
+ c->SetFillColor(0);
+ c->SetBorderMode(0);
+ c->SetBorderSize(0);
+
+
+ TArrayF measurements(6);
+ std::cout << "Measurements are: " << std::flush;
+ for (Int_t i = 0; i < measurements.fN; i++) {
+ measurements[i] = gRandom->Uniform(0,1);
+ std::cout << measurements[i] << " " << std::flush;
+ }
+ std::cout << std::endl;
+
+ Float_t last = 0;
+ Float_t B = 5;
+ TH2* frame = new TH2F("frame", "Frame", measurements.fN, 0,
+ measurements.fN, 10, 0, 1.1);
+ frame->Draw();
+ for (Int_t i = 0; i < measurements.fN; i++) {
+ TF1* f = new TF1("f", "[2] + exp(-[1] * (x - [3])) * ([0] - [2])",
+ i, i + 1);
+ f->SetParameter(3, i);
+ f->SetParameter(1, B);
+ f->SetParameter(2, measurements[i]);
+ f->SetParameter(0, last);
+
+ if (measurements[i] > last) {
+ // f = new TF1("f", "[0] * (1 - exp(-[1] * (x - [2]))) + [3]", i, i+1);
+ // f->SetParameters(measurements[i] - last, B, i, last);
+ // f->SetParameter(0, measurements[i]);
+ // f->SetParameter(2, last);
+ }
+ else {
+ // f->SetParameter(2, measurements[i]);
+ // f->SetParameter(0, last);
+ // f = new TF1("f", "[0] * (exp(-[1] * (x - [2]))) + [3]", i, i+1);
+ // f->SetParameters(last - measurements[i], B, i, measurements[i]);
+ }
+ f->Draw("same");
+ last = measurements[i];
+ }
+}
+
+
+
--- /dev/null
+void ViewPIPE()
+{
+ gMC->Gsatt("QBPM","seen",0);
+ // gMC->Gsatt("QB11","seen",0);
+ // gMC->Gsatt("QB00","seen",1);
+ // gMC->Gsatt("QB02","seen",1);
+ // gMC->Gsatt("QB01","seen",1);
+ // gMC->Gsatt("QB04","seen",1);
+ // gMC->Gsatt("QB03","seen",1);
+ // gMC->Gsatt("QB05","seen",1);
+ // gMC->Gsatt("QB06","seen",1);
+ // gMC->Gsatt("QB08","seen",1);
+ // gMC->Gsatt("QB10","seen",1);
+ // gMC->Gsatt("QB07","seen",1);
+ // gMC->Gsatt("QB30","seen",1);
+ // gMC->Gsatt("QB27","seen",1);
+ // gMC->Gsatt("QB26","seen",1);
+ // gMC->Gsatt("QB25","seen",1);
+ // gMC->Gsatt("QB20","seen",1);
+ // gMC->Gsatt("QB23","seen",1);
+ // gMC->Gsatt("QB22","seen",1);
+ // gMC->Gsatt("QB19","seen",1);
+ // gMC->Gsatt("QB18","seen",1);
+ // gMC->Gsatt("QB21","seen",1);
+ // gMC->Gsatt("QB15","seen",1);
+ // gMC->Gsatt("QB16","seen",1);
+ // gMC->Gsatt("QB17","seen",1);
+ // gMC->Gsatt("QB14","seen",1);
+ // gMC->Gsatt("QB13","seen",1);
+ // gMC->Gsatt("QB12","seen",1);
+ // gMC->Gsatt("QB31","seen",1);
+ // gMC->Gsatt("QB32","seen",1);
+ // gMC->Gsatt("QIPM","seen",0);
+
+ gMC->Gsatt("QI32","seen",1);
+ gMC->Gsatt("QI33","seen",1);
+ gMC->Gsatt("QI34","seen",1);
+ gMC->Gsatt("QI35","seen",1);
+ gMC->Gsatt("QI42","seen",1);
+ gMC->Gsatt("QI43","seen",1);
+ gMC->Gsatt("QI42","seen",1);
+
+ gMC->Gsatt("QFA0","seen",1);
+
+ gMC->Gsatt("QB24","seen",1);
+ gMC->Gsatt("QA24","seen",1);
+ gMC->Gsatt("QB28","seen",1);
+
+ gMC->Gsatt("QBE0","seen",0);
+ gMC->Gsatt("QBEP","seen",1);
+ gMC->Gsatt("QBEM","seen",1);
+ gMC->Gsatt("QBEW","seen",1);
+ gMC->Gsatt("QBEU","seen",1);
+
+ gMC->Gsatt("QB29","seen",0);
+ gMC->Gsatt("QP29","seen",1);
+ gMC->Gsatt("QS29","seen",1);
+ gMC->Gsatt("QF29","seen",1);
+
+ gMC->Gsatt("QBAB","seen",1);
+ gMC->Gsatt("QBBE","seen",1);
+
+ gMC->Gsatt("QBSR","seen",1);
+ gMC->Gsatt("QBSS","seen",1);
+
+ gMC->Gsatt("QBVA","seen",1);
+}
--- /dev/null
+//
+// Small script that I used to make some intial testing of the wafer
+// layout and geometry.
+//
+// Christian
+//
+TObjArray*
+WaferParameters(const char c)
+{
+ double dl;
+ double dh;
+ double r = 134 / 2;
+ double theta;
+ switch (c) {
+ case 'i':
+ dl = 43;
+ dh = 172;
+ theta = 18;
+ break;
+ case 'o':
+ dl = 156;
+ dh = 280;
+ theta = 9;
+ break;
+ default:
+ cerr << "Unknown wafer type: " << c << endl;
+ return;
+ }
+
+
+ double tan_theta = tan(theta * TMath::Pi() / 180.);
+ double tan_theta2 = pow(tan_theta,2);
+ double r2 = pow(r,2);
+ double y_A = tan_theta * dl;
+ double x_D = dl + sqrt(r2 - tan_theta2 * pow(dl,2));
+ double x_D_2 = dl - sqrt(r2 - tan_theta2 * pow(dl,2));
+
+ double y_B = sqrt(r2 - pow(dh,2) + 2 * dh * x_D - pow(x_D,2));
+ double x_C = (x_D + sqrt(-tan_theta2 * pow(x_D,2) + r2
+ + r2 * tan_theta2)) / (1 + tan_theta2);
+ double y_C = tan_theta * x_C;
+
+ cout << "A: (" << dl << "," << y_A << ")" << endl;
+ cout << "B: (" << dh << "," << y_B << ")" << endl;
+ cout << "C: (" << x_C << "," << y_C << ")" << endl;
+ cout << "D: (" << x_D << ",0)" << endl;
+
+ cout << "Recentred at D:" << endl;
+ cout << "A: (" << dl - x_D << "," << y_A << ")" << endl;
+ cout << "B: (" << dh - x_D << "," << y_B << ")" << endl;
+ cout << "C: (" << x_C - x_D << "," << y_C << ")" << endl;
+
+ TObjArray* verticies = new TObjArray(6);
+ verticies->AddAt(new TVector2(dl, y_A), 0);
+ verticies->AddAt(new TVector2(x_C, y_C), 1);
+ verticies->AddAt(new TVector2(dh, y_B), 2);
+ verticies->AddAt(new TVector2(dh, -y_B), 3);
+ verticies->AddAt(new TVector2(x_C, -y_C), 4);
+ verticies->AddAt(new TVector2(dl, -y_A), 5);
+
+ return verticies;
+}
+
+void
+Wafer()
+{
+ TCanvas* can = new TCanvas("can", "c", 400, 600);
+ can->SetBorderMode(0);
+ can->SetFillColor(0);
+
+ TGeometry* g = new TGeometry("g", "G");
+ TShape* topShape = new TBRIK("top", "top", "", 100, 100, 100);
+ TNode* topNode = new TNode("top", "top", "top", 0, 0, 0);
+ topNode->SetLineWidth(0);
+ topNode->SetVisibility(0);
+
+ TShape* virtualShape = new TTUBS("virtual", "Virtual", "",
+ 43, 172, 1, -18, 18);
+
+ TObjArray* v = WaferParameters('i');
+ TXTRU* moduleShape = new TXTRU("module", "module", "", 6, 2);
+ for (Int_t i = 0; i < 6; i++) {
+ TVector2* vv = static_cast<TVector2*>(v->At(i));
+ moduleShape->DefineVertex(i, vv->X(), vv->Y());
+ }
+ moduleShape->DefineSection(0, -1, 1, 0, 0);
+ moduleShape->DefineSection(1, 1, 1, 0, 0);
+
+ for (Int_t i = 0; i < 10; i++) {
+ topNode->cd();
+ Double_t theta = 36 * i;
+ Double_t z = (i % 2 ? +5 : -5);
+ TRotMatrix* rot = new TRotMatrix(Form("rotation%02d", i), "Rotation",
+ 90, theta, 90,
+ fmod(90 + theta, 360), 0, 0);
+ TNode* moduleNode = new TNode(Form("module%02d", i),
+ "Module", moduleShape, 0, 0, z,
+ rot);
+ if (i == 0) {
+ moduleNode->SetFillColor(2);
+ moduleNode->SetLineColor(2);
+ moduleNode->SetLineWidth(2);
+ }
+ }
+ g->Draw();
+ TView* view = can->GetView();
+ view->SetPerspective();
+ Int_t a;
+ view->SetView(1.81208, 66.6725, 90, a);
+ view->Zoom();
+ view->Zoom();
+ view->Zoom();
+ can->Modified();
+ can->cd();
+
+ can->Print("fmd_module1.gif");
+ // std::cout << "Waiting ..." << std::flush;
+ // Char_t c = std::cin.get();
+
+ topNode->cd();
+ TNode* virtualNode = new TNode("virtual", "Virtual",
+ virtualShape, 0, 0, -5);
+ virtualNode->SetLineColor(3);
+ virtualNode->SetLineWidth(2);
+ virtualNode->SetLineStyle(2);
+ g->Draw();
+ view->SetPerspective();
+ view->SetView(1.81208, 66.6725, 90, a);
+ view->Zoom();
+ view->Zoom();
+ view->Zoom();
+ can->Modified();
+ can->cd();
+ can->Print("fmd_module2.gif");
+
+}
+
+
+
+
+
+
--- /dev/null
+#!/bin/bash
+#
+# Shell script to do all the LEGO plots
+#
+for i in Inner ITS PIPE FMD Nothing ; do
+ aliroot -l -b -q FMD/scripts/MakeLego.C\(\"$i\"\)
+done
+
+root -l -q FMD/scripts/DrawLego.C
--- /dev/null
+#!/bin/bash
+CURDIR=`pwd`
+cd $ALICE_ROOT
+echo 'Making sure that TFluka is up to date...'
+make all-TFluka
+cd $CURDIR
+
+# Make working directory
+rm -rf fluka
+mkdir -p fluka
+cd fluka
+
+# Make a peg directory
+mkdir -p peg
+
+# Link here some special Fluka files needed
+ln -s $FLUPRO/xnloan.dat .
+ln -s $FLUPRO/sigmapi.bin .
+ln -s $FLUPRO/nuclear.bin .
+ln -s $FLUPRO/neuxsc_72.bin neuxsc.bin
+ln -s $FLUPRO/fluodt.dat .
+ln -s $FLUPRO/elasct.bin .
+
+# Copy the random seed
+cp $FLUPRO/random.dat old.seed
+
+# Give some meaningfull name to the output
+ln -s fluka.out fort.11
+
+# Link the pemf and input file for alice
+# This is wrong:
+# ln -s $ALICE_ROOT/TFluka/input/FlukaVmc.pemf .
+# Maybe
+cp $ALICE_ROOT/TFluka/input/alice.pemf FlukaVmc.pemf
+
+#Link FlukaConfig.C as Config.C
+cp $ALICE_ROOT/FMD/scripts/ConfigInner.C .
+cp $ALICE_ROOT/.rootrc .
+# echo 'Execute: gAlice->Init() OR gAlice->RunMC() at the ROOT prompt'
+# Launch aliroot
+aliroot -l # -b -q ../runIt.C > run.log 2>&1
+
+cd $CURDIR
--- /dev/null
+#!/bin/bash
+CURDIR=`pwd`
+cd $ALICE_ROOT
+echo 'Making sure that AliROOT is up to date...'
+make
+cd $CURDIR
+
+# Make working directory
+rm -rf geant321
+mkdir -p geant321
+cd geant321
+
+
+#Link FlukaConfig.C as Config.C
+cp $ALICE_ROOT/FMD/Config.C .
+cp $ALICE_ROOT/.rootrc .
+# echo 'Execute: gAlice->Init() OR gAlice->RunMC() at the ROOT prompt'
+# Launch aliroot
+aliroot -l -b -q ../runIt.C > run.log 2>&1
+
+cd $CURDIR
git://git.uio.no / u/mrichter/AliRoot.git / commitdiff