--- /dev/null
+/**************************************************************************
+ * Copyright(c) 2004, 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$ */
+/** @file AliFMDDigit.cxx
+ @author Christian Holm Christensen <cholm@nbi.dk>
+ @date Mon Mar 27 12:37:41 2006
+ @brief Digits for the FMD
+*/
+//////////////////////////////////////////////////////////////////////
+//
+// Class that holds an FMD index. That is, it holds the detector
+// coordinates for a given strip:
+//
+// Variable | Type | Range | Description
+// ---------+----------+---------+------------------
+// detector | UShort_t | 1-3 | Detector number
+// ring | Char_t | 'I'/'O' | Ring identifier
+// sector | UShort_t | 0-39 | Sector number
+// strip | UShort_t | 0-511 | Strip number
+//
+//////////////////////////////////////////////////////////////////////
+
+#include "AliFMDIndex.h" // ALIFMDINDEX_H
+#include "Riostream.h" // ROOT_Riostream
+#include <TString.h> // ROOT_TString
+#include <AliFMDMap.h>
+
+//====================================================================
+ClassImp(AliFMDIndex)
+#if 0
+ ; // This is here to keep Emacs from indenting the next line
+#endif
+
+//____________________________________________________________________
+AliFMDIndex::AliFMDIndex()
+ : fDetector(0),
+ fRing('\0'),
+ fSector(0),
+ fStrip(0),
+ fHash(-1)
+{}
+
+//____________________________________________________________________
+AliFMDIndex::AliFMDIndex(const AliFMDIndex& o)
+ : fDetector(o.fDetector),
+ fRing(o.fRing),
+ fSector(o.fSector),
+ fStrip(o.fStrip),
+ fHash(o.fHash)
+{
+ // Copy constructor
+}
+
+//____________________________________________________________________
+AliFMDIndex::AliFMDIndex(UShort_t detector,
+ Char_t ring,
+ UShort_t sector,
+ UShort_t strip)
+ : fDetector(detector),
+ fRing(ring),
+ fSector(sector),
+ fStrip(strip),
+ fHash(-1)
+{
+ //
+ // Creates a base data digit object
+ //
+ // Parameters
+ //
+ // detector Detector # (1, 2, or 3)
+ // ring Ring ID ('I' or 'O')
+ // sector Sector # (For inner/outer rings: 0-19/0-39)
+ // strip Strip # (For inner/outer rings: 0-511/0-255)
+}
+
+//____________________________________________________________________
+AliFMDIndex&
+AliFMDIndex::operator=(const AliFMDIndex& o)
+{
+ // Assignment operator
+ fDetector = o.fDetector;
+ fRing = o.fRing;
+ fSector = o.fSector;
+ fStrip = o.fStrip;
+ fHash = o.fHash;
+ return *this;
+}
+
+//____________________________________________________________________
+Int_t
+AliFMDIndex::Hash() const
+{
+ if (fHash < 0) {
+ size_t ringi = (fRing == 'I' || fRing == 'i' ? 0 : 1);
+ fHash = (fStrip +
+ AliFMDMap::kMaxStrips *
+ (fSector + AliFMDMap::kMaxSectors *
+ (ringi + AliFMDMap::kMaxRings * (fDetector-1))));
+ }
+ return fHash;
+}
+
+
+//____________________________________________________________________
+void
+AliFMDIndex::Print(Option_t* /* option*/) const
+{
+ // Print digit to standard out
+ cout << Name() << flush;
+}
+
+//____________________________________________________________________
+const char*
+AliFMDIndex::Name() const
+{
+ if (fName.IsNull())
+ fName = Form("FMD%d%c[%2d,%3d]", fDetector, fRing, fSector, fStrip);
+ return fName.Data();
+}
+
+//====================================================================
+ClassImp(AliFMDObjIndex)
+#if 0
+ ; // This is here to keep Emacs from indenting the next line
+#endif
+
+//____________________________________________________________________
+Int_t
+AliFMDObjIndex::Compare(const TObject* o) const
+{
+ const AliFMDObjIndex* a = dynamic_cast<const AliFMDObjIndex*>(o);
+ if (!a) {
+ Fatal("Compare",
+ "trying to compare to something not a AliFMDObjIndex object, "
+ "but a %s object", o->ClassName());
+ return 0;
+ }
+ if (this->operator<(*a)) return -1;
+ if (this->operator==(*a)) return 0;
+ return 1;
+}
+
+//____________________________________________________________________
+//
+// EOF
+//
--- /dev/null
+#ifndef ALIFMDINDEX_H
+#define ALIFMDINDEX_H
+/** @file AliFMDIndex.h
+ @author Christian Holm Christensen <cholm@nbi.dk>
+ @date Mon Mar 27 12:37:41 2006
+ @brief FMD detector coordinates
+*/
+//___________________________________________________________________
+//
+// Class that holds an FMD index. That is, it holds the detector
+// coordinates for a given strip:
+//
+// Variable | Type | Range | Description
+// ---------+----------+---------+------------------
+// detector | UShort_t | 1-3 | Detector number
+// ring | Char_t | 'I'/'O' | Ring identifier
+// sector | UShort_t | 0-39 | Sector number
+// strip | UShort_t | 0-511 | Strip number
+//
+#ifndef ROOT_Rtypes
+# include <Rtypes.h>
+#endif
+#ifndef ROOT_TObject
+# include <TObject.h>
+#endif
+#ifndef ROOT_TString
+# include <TString.h>
+#endif
+#include <iosfwd>
+
+//____________________________________________________________________
+/** @class AliFMDIndex AliFMDIndex.h <FMD/AliFMDIndex.h>
+ @brief FMD detector coordinates
+ @ingroup FMD_base
+ */
+class AliFMDIndex
+{
+public:
+ /** CTOR */
+ AliFMDIndex();
+ /** Copy CTOR
+ @param o Object to copy from */
+ AliFMDIndex(const AliFMDIndex& o);
+ /** Constrctor
+ @param detector Detector
+ @param ring Ring
+ @param sector Sector
+ @param strip Strip */
+ AliFMDIndex(UShort_t detector,
+ Char_t ring='\0',
+ UShort_t sector=0,
+ UShort_t strip=0);
+ /** Assignment operator
+ @param o Object to assign from
+ @return Reference to this object */
+ AliFMDIndex& operator=(const AliFMDIndex& o);
+ /** Comparison operator
+ @param o Object to compare to
+ @return @c true if these refer to the same index */
+ bool operator==(const AliFMDIndex& o) const;
+ /** Comparison operator
+ @param o Object to compare to
+ @return @c true if this is smaller than @a o */
+ bool operator<(const AliFMDIndex& o) const;
+ /** DTOR */
+ virtual ~AliFMDIndex() {}
+ /** @return Detector # */
+ UShort_t Detector() const { return fDetector; }
+ /** @return Ring ID */
+ Char_t Ring() const { return fRing; }
+ /** @return sector # */
+ UShort_t Sector() const { return fSector; }
+ /** @return strip # */
+ UShort_t Strip() const { return fStrip; }
+ /** @param x Detector # */
+ void SetDetector(UShort_t x) { fHash = -1; fDetector = x; }
+ /** @param x Ring ID */
+ void SetRing(Char_t x) { fHash = -1; fRing = x; }
+ /** @param x sector # */
+ void SetSector(UShort_t x) { fHash = -1; fSector = x; }
+ /** @param x strip # */
+ void SetStrip(UShort_t x) { fHash = -1; fStrip = x; }
+ /** Print information
+ @param opt Not used */
+ virtual void Print(Option_t* opt="") const;
+ /** @return Name */
+ const char* Name() const;
+protected:
+ Int_t Hash() const;
+ UShort_t fDetector; // (Sub) Detector # (1,2, or 3)
+ Char_t fRing; // Ring ID ('I' or 'O')
+ UShort_t fSector; // Sector # (phi division)
+ UShort_t fStrip; // Strip # (radial division)
+ mutable TString fName; //! Cached name
+ mutable Int_t fHash; //! Cached hash value
+ ClassDef(AliFMDIndex, 1) // Base class for FMD digits
+};
+
+//____________________________________________________________________
+class AliFMDObjIndex : public TObject, public AliFMDIndex
+{
+public:
+ /** CTOR */
+ AliFMDObjIndex() {}
+ /** Copy CTOR
+ @param o Object to copy from */
+ AliFMDObjIndex(const AliFMDObjIndex& o) : TObject(o), AliFMDIndex(o) {}
+ /** Construct from a pure index
+ @param o Object to copy from */
+ explicit AliFMDObjIndex(const AliFMDIndex& o) : AliFMDIndex(o) {}
+ /** Constrctor
+ @param detector Detector
+ @param ring Ring
+ @param sector Sector
+ @param strip Strip */
+ AliFMDObjIndex(UShort_t detector,
+ Char_t ring='\0',
+ UShort_t sector=0,
+ UShort_t strip=0)
+ : AliFMDIndex(detector, ring, sector, strip)
+ {}
+ /** DTOR */
+ virtual ~AliFMDObjIndex() {}
+ /** @return name */
+ virtual const char* GetName() const { return AliFMDIndex::Name(); }
+ /** sort compare for TCollection's
+ @param o Object to compare to
+ @return -1 if this is @e smaller than @a o, 0 if @e equal to
+ @a o, and 1 if this is @e larger than @a o */
+ virtual Int_t Compare(const TObject* o) const;
+ /** @return always true */
+ Bool_t IsSortable() const { return kTRUE; }
+ ClassDef(AliFMDObjIndex, 1) // Base class for FMD digits
+};
+
+//____________________________________________________________________
+inline
+bool
+AliFMDIndex::operator==(const AliFMDIndex& o) const
+{
+ return (o.Hash() == Hash());
+}
+
+//____________________________________________________________________
+inline
+bool
+AliFMDIndex::operator<(const AliFMDIndex& rhs) const
+{
+ return (Hash() < rhs.Hash());
+}
+
+#if 0
+//____________________________________________________________________
+inline
+bool
+operator<(const AliFMDIndex& lhs, const AliFMDIndex& rhs)
+{
+ return (lhs.Detector() < rhs.Detector() ? true :
+ (lhs.Ring() < rhs.Ring() ? true :
+ (lhs.Sector() < rhs.Sector() ? true :
+ (lhs.Strip() < rhs.Strip() ? true : false))));
+}
+#endif
+#endif
+//____________________________________________________________________
+//
+// Local Variables:
+// mode: C++
+// End:
+//
+//
+// EOF
+//
// Get sample rate
AliFMDParameters* pars = AliFMDParameters::Instance();
AliFMDRawStream input(fReader);
- // Select FMD DDL's
- fReader->Select(AliFMDParameters::kBaseDDL>>8);
UShort_t stripMin = 0;
UShort_t stripMax = 127;
: AliAltroRawStream(reader)
{
fNoAltroMapping = kFALSE;
+ // Select FMD DDL's
+ SelectRawData(AliFMDParameters::kBaseDDL>>8);
}
//_____________________________________________________________________________
AliFMDRawStream::ReadChannel(UInt_t& ddl, UInt_t& addr,
UInt_t& len, UShort_t* data)
{
- UInt_t prevddl = 0;
Int_t l = 0;
static Int_t last = 0xFFFF; // 0xFFFF means signal is used
Bool_t next = kTRUE;
AliDebug(15, Form("New hardware address, was 0x%x, now 0x%x",
GetPrevHWAddress(), GetHWAddress()));
addr = GetPrevHWAddress();
- ddl = AliFMDParameters::kBaseDDL + prevddl;
- len = l+1;
+ ddl = AliFMDParameters::kBaseDDL + GetPrevDDLNumber();
+ len = l+1; // Need to add one - l points to last valid index
last = signal;
break;
}
}
- prevddl = GetPrevDDLNumber();
Int_t t = GetTime();
l = TMath::Max(l, t);
data[t] = signal;
# Configuration options related to the dot tool
#---------------------------------------------------------------------------
CLASS_DIAGRAMS = YES
-HIDE_UNDOC_RELATIONS = YES
+HIDE_UNDOC_RELATIONS = NO
HAVE_DOT = YES
CLASS_GRAPH = YES
COLLABORATION_GRAPH = YES
-GROUP_GRAPHS = YES
+GROUP_GRAPHS = NO
UML_LOOK = YES
TEMPLATE_RELATIONS = NO
-INCLUDE_GRAPH = YES
-INCLUDED_BY_GRAPH = YES
+INCLUDE_GRAPH = NO
+INCLUDED_BY_GRAPH = NO
CALL_GRAPH = NO
GRAPHICAL_HIERARCHY = YES
-DIRECTORY_GRAPH = YES
+DIRECTORY_GRAPH = NO
DOT_IMAGE_FORMAT = png
DOT_PATH =
DOTFILE_DIRS =
#pragma link off all classes;
#pragma link off all functions;
+#pragma link C++ class AliFMDIndex+;
+#pragma link C++ class AliFMDObjIndex+;
#pragma link C++ class AliFMDBaseDigit+;
#pragma link C++ class AliFMDDigit+;
#pragma link C++ class AliFMDSDigit+;
cdb->SetDefaultStorage("local://$ALICE_ROOT");
AliLog::SetModuleDebugLevel("FMD", 2);
AliReconstruction rec;
+ AliCDBEntry* align = cdb->Get("FMD/Align/Data");
+ if (align) {
+ TClonesArray* array = dynamic_cast<TClonesArray*>(align->GetObject());
+ if (array) {
+ Int_t nAlign = array->GetEntries();
+ for (Int_t i = 0; i < nAlign; i++) {
+ AliAlignObjAngles* a = static_cast<AliAlignObjAngles*>(array->At(i));
+ if (!a->ApplyToGeometry()) {
+ Warning("ApplyAlignement", "Failed to apply alignment to %s",
+ a->GetVolPath());
+ }
+ }
+ }
+ }
rec.SetRunLocalReconstruction("FMD");
rec.SetRunVertexFinder(kFALSE);
rec.SetRunTracking("");
#
# $Id$
-SRCS = AliFMDDigit.cxx \
+SRCS = AliFMDIndex.cxx \
+ AliFMDDigit.cxx \
AliFMDBoolMap.cxx \
AliFMDUShortMap.cxx \
AliFMDCalibPedestal.cxx \
--- /dev/null
+//____________________________________________________________________
+//
+//
+// $Id$
+//
+// Script I used for rapid prototyping of the FMD3 geometry - in
+// particular the support cone
+//
+/** @defgroup geo_geom Simple geometry
+ @ingroup FMD_script
+*/
+#ifndef __CINT___
+# ifndef ROOT_TGeoVolume
+# include <TGeoVolume.h>
+# endif
+# ifndef ROOT_TGeoMaterial
+# include <TGeoMaterial.h>
+# endif
+# ifndef ROOT_TGeoMedium
+# include <TGeoMedium.h>
+# endif
+# ifndef ROOT_TGeoMatrix
+# include <TGeoMatrix.h>
+# endif
+# ifndef ROOT_TGeoXtru
+# include <TGeoXtru.h>
+# endif
+# ifndef ROOT_TGeoPcon
+# include <TGeoPcon.h>
+# endif
+# ifndef ROOT_TGeoTube
+# include <TGeoTube.h>
+# endif
+# ifndef ROOT_TGeoManager
+# include <TGeoManager.h>
+# endif
+# ifndef ROOT_TGeoPhysicalNode
+# include <TGeoPhysicalNode.h>
+# endif
+# ifndef ROOT_TVector3
+# include <TVector3.h>
+# endif
+# ifndef ROOT_TString
+# include <TString.h>
+# endif
+# ifndef ROOT_TRandom
+# include <TRandom.h>
+# endif
+# ifndef ROOT_TBrowser
+# include <TBrowser.h>
+# endif
+# ifndef ROOT_TCanvas
+# include <TCanvas.h>
+# endif
+# ifndef __CSTDARG__
+# include <cstdarg>
+# endif
+# ifndef __IOSTREAM__
+# include <iostream>
+# endif
+#endif
+
+//____________________________________________________________________
+/** @brief Simple geometry
+ @ingroup geo_geom
+ @code
+ gSystem->Load("libPhysics.so");
+ gSystem->Load("libGeom.so");
+ gROOT->LoadMacro("GeoGeometry.C+");
+ Geometry g;
+ g.Exec();
+ @endcode
+ */
+class Geometry
+{
+public:
+ /** Constructor */
+ Geometry();
+ /** Destructor */
+ virtual ~Geometry() {}
+ /** Initialize */
+ virtual void Initialize();
+ /** Register */
+ virtual void Register();
+ /** @e Do-It member function */
+ virtual void Exec();
+ /** Try to align */
+ virtual void Align();
+ /** Convert detector coordinates to spatial coordinates
+ @param sector Sector number
+ @param strip Strip number
+ @param xyz Spatial coordinates
+ */
+ virtual void Detector2XYZ(UInt_t sector, UInt_t strip, TVector3& xyz);
+ /** Debug level */
+ static Int_t fgDebug;
+ /** Print debugging messages
+ @param lvl Acceptance level
+ @param where Where it happened
+ @param format Message format. */
+ static void Debug(Int_t lvl, const char* where, const char* format, ...)
+ {
+ if (lvl > fgDebug) return;
+ static char buf[2048];
+ va_list ap;
+ va_start(ap, format);
+ vsnprintf(buf, 2048, format, ap);
+ std::cout << "D: " << where << ": " << buf << std::endl;
+ va_end(ap);
+ }
+protected:
+ /** List of translations */
+ TObjArray* fMatricies;
+ /** Shape parameter */
+ TVector2 fA;
+ /** Shape parameter */
+ TVector2 fB;
+ /** Shape parameter */
+ TVector2 fC;
+ /** Spacing between sensor and hybrid */
+ Double_t fSpacer;
+ /** Thickness of aluminum plates in honey comb */
+ Double_t fAlThickness;
+ /** Width of bonding pads */
+ Double_t fBondingWidth;
+ /** chip layer thickenss */
+ Double_t fChipThickness;
+ /** Copper layer thickness */
+ Double_t fCopperThickness;
+ /** Upper radious */
+ Double_t fHighR;
+ /** Thickness of honey comb */
+ Double_t fHoneycombThickness;
+ /** Inner honey comb inner radius */
+ Double_t fInnerHoneyHighR;
+ /** Inner honey comb outer radius */
+ Double_t fInnerHoneyLowR;
+ /** Z coordinate of inner ring */
+ Double_t fInnerZ;
+ /** Length of support legs */
+ Double_t fLegLength;
+ /** Offset from edge of legs */
+ Double_t fLegOffset;
+ /** Radius of support legs */
+ Double_t fLegRadius;
+ /** Inner radius */
+ Double_t fLowR;
+ /** Spacing between modules */
+ Double_t fModuleSpacing;
+ /** Number of strps */
+ Int_t fNStrips;
+ /** Outer honey comb inner radius */
+ Double_t fOuterHoneyHighR;
+ /** Outer honey comb outer radius */
+ Double_t fOuterHoneyLowR;
+ /** Z coordinate of outer */
+ Double_t fOuterZ;
+ /** Thickness of print board */
+ Double_t fPrintboardThickness;
+ /** Cache of ring depth */
+ Double_t fRingDepth;
+ /** Thickness of silicon sensor */
+ Double_t fSiThickness;
+ /** ??? */
+ Double_t fSpacerHeight;
+ /** Opening angle of sensor */
+ Double_t fTheta;
+ /** Radius of wafer the sensor is cut out of */
+ Double_t fWaferRadius;
+ /** Air tracking medium */
+ TGeoMedium* fAir;
+ /** Aluminum tracking medium */
+ TGeoMedium* fAl;
+ /** Carbon tracking medium */
+ TGeoMedium* fCarbon;
+ /** Chip tracking medium */
+ TGeoMedium* fChip;
+ /** Copper tracking medium */
+ TGeoMedium* fCopper;
+ /** Kapton tracking medium */
+ TGeoMedium* fKapton;
+ /** PCB tracking medium */
+ TGeoMedium* fPCB;
+ /** Plastic tracking medium */
+ TGeoMedium* fPlastic;
+ /** Active silicon tracking medium */
+ TGeoMedium* fSi;
+ /** Vacuum tracking medium */
+ TGeoMedium* fVacuum;
+ /** Use assemblies */
+ Bool_t fDoubleAssembly;
+ /** Make a detailed geometry */
+ void MakeDetailed(TGeoVolume* mother);
+};
+
+//____________________________________________________________________
+Int_t Geometry::fgDebug = 10;
+
+//____________________________________________________________________
+Geometry::Geometry()
+ : fMatricies(0),
+ fA( 4.3000, 1.3972),
+ fB(17.2000, 2.1452),
+ fC(15.3327, 4.9819)
+{
+ fBondingWidth = 0.5000;
+ fWaferRadius = 6.7000;
+ fSiThickness = 0.0300;
+ fLowR = 4.3000;
+ fHighR = 17.2000;
+ fTheta = 18.0000;
+ fNStrips = 512;
+ fRingDepth = 2.1600;
+ fLegRadius = 0.5000;
+ fLegLength = 1.0000;
+ fLegOffset = 2.0000;
+ fModuleSpacing = 1.0000;
+ fPrintboardThickness = 0.1000;
+ fSpacerHeight = 0.0300;
+
+ fInnerZ = 83.4000;
+ fOuterZ = 75.2000;
+ fHoneycombThickness = 1.0000;
+ fAlThickness = 0.1000;
+ fInnerHoneyLowR = 5.3000;
+ fInnerHoneyHighR = 18.2000;
+ fOuterHoneyLowR = 16.6000;
+ fOuterHoneyHighR = 29.0000;
+
+ fCopperThickness = 0.01;
+ fChipThickness = 0.01;
+ fAlThickness = 0.10;
+ fHoneycombThickness = 1.00;
+ fSpacer = 0.10;
+ fLegLength = 1.00;
+}
+
+//____________________________________________________________________
+void
+Geometry::Exec()
+{
+ Initialize();
+ Register();
+ // gGeoManager->DefaultColors();
+ gGeoManager->ViewLeaves(kFALSE);
+ TCanvas* c = new TCanvas;
+ c->SetFillColor(0);
+ gGeoManager->GetTopVolume()->Draw();
+ new TBrowser("b", "Browser");
+}
+
+//____________________________________________________________________
+void
+Geometry::Initialize()
+{
+ Double_t mMax = 10;
+ Double_t mType = 2;
+
+ // Air
+ Double_t pAir[] = { 0., mType, mMax, 1., .001, 1., .001, .001 };
+ TGeoMixture* fmdAir = new TGeoMixture("FMD air", 4, .00120479);
+ fmdAir->DefineElement(0, 12.0107, 6., 0.000124);
+ fmdAir->DefineElement(1, 14.0067, 7., 0.755267);
+ fmdAir->DefineElement(2, 15.9994, 8., 0.231781);
+ fmdAir->DefineElement(3, 39.948, 18., 0.012827);
+ fmdAir->SetTransparency('0');
+ fAir = new TGeoMedium("FMD Air", 1, fmdAir, pAir);
+
+ // Silicon
+ Double_t pSi[] = { 1., mType, mMax, 1., .001, 1., .001, .001 };
+ TGeoMaterial* fmdSi = new TGeoMaterial("FMD Si", 28.0855, 14, 2.33);
+ fmdSi->SetFillColor(2);
+ fSi = new TGeoMedium("FMD Si", 1, fmdSi, pSi);
+
+ // Vacumm
+ Double_t pVacuum[] = { 0., mType, mMax, 10, .01, .1, .003, .003 };
+ TGeoMaterial* fmdVacuum = new TGeoMaterial("FMD Vacuum",1e-16,1e-16,1e-16);
+ fmdVacuum->SetTransparency('0');
+ fVacuum = new TGeoMedium("FMD Vacuum", 1, fmdVacuum,pVacuum);
+
+
+ // PCB
+ Double_t pPCB[] = { 0., mType, mMax, 1., .001, 1., .001, .001 };
+ TGeoMixture* fmdPCB = new TGeoMixture("FMD PCB", 14, 1.8);
+ fmdPCB->DefineElement(0, 28.0855, 14, 0.15144894);
+ fmdPCB->DefineElement(1, 40.078, 20, 0.08147477);
+ fmdPCB->DefineElement(2, 26.981538, 13, 0.04128158);
+ fmdPCB->DefineElement(3, 24.305, 12, 0.00904554);
+ fmdPCB->DefineElement(4, 10.811, 5, 0.01397570);
+ fmdPCB->DefineElement(5, 47.867, 22, 0.00287685);
+ fmdPCB->DefineElement(6, 22.98977, 11, 0.00445114);
+ fmdPCB->DefineElement(7, 39.0983, 19, 0.00498089);
+ fmdPCB->DefineElement(8, 55.845, 26, 0.00209828);
+ fmdPCB->DefineElement(9, 18.9984, 9, 0.00420000);
+ fmdPCB->DefineElement(10, 15.9994, 8, 0.36043788);
+ fmdPCB->DefineElement(11, 12.0107, 6, 0.27529425);
+ fmdPCB->DefineElement(12, 14.0067, 7, 0.01415852);
+ fmdPCB->DefineElement(13, 1.00794, 1, 0.03427566);
+ fmdPCB->SetFillColor(7);
+ fPCB = new TGeoMedium("FMD PCB", 1, fmdPCB, pPCB);
+
+ // Chip
+ Double_t pChip[] = { 0., mType, mMax, 10., .01, 1., .003, .003 };
+ TGeoMixture* fmdChip = new TGeoMixture("FMD Chip", 6, 2.36436);
+ fmdChip->DefineElement(0, 12.0107, 6, 0.039730642);
+ fmdChip->DefineElement(1, 14.0067, 7, 0.001396798);
+ fmdChip->DefineElement(2, 15.9994, 8, 0.01169634);
+ fmdChip->DefineElement(3, 1.00794, 1, 0.004367771);
+ fmdChip->DefineElement(4, 28.0855, 14, 0.844665);
+ fmdChip->DefineElement(5, 107.8682, 47, 0.0981434490);
+ fmdChip->SetFillColor(4);
+ fChip = new TGeoMedium("FMD Chip", 1, fmdChip, pChip);
+
+ // Carbon
+ Double_t pC[] = { 0., mType, mMax, 10., .01, 1., .003, .003 };
+ TGeoMaterial* fmdC = new TGeoMaterial("FMD C", 12.011, 6, 2.265);
+ fmdC->SetFillColor(6);
+ fCarbon = new TGeoMedium("FMD C", 1, fmdC, pC);
+
+ // Kapton (inside of Honeycombs)
+ Double_t pKapton[] = { 0., mType, mMax, 1., .001, 1., .001, .001 };
+ TGeoMaterial* fmdKapton = new TGeoMaterial("FMD Kapton", 12.011, 6.,0.01);
+ fmdKapton->SetFillColor(3);
+ fKapton = new TGeoMedium("FMD Kapton",1,fmdKapton,pKapton);
+
+ // Plastic
+ Double_t pPlastic[] = { 0., mType, mMax, 10., .01, 1., .003, .003 };
+ TGeoMixture* fmdPlastic = new TGeoMixture("FMD Plastic", 2, 1.03);
+ fmdPlastic->DefineElement(0, 1.01, 1, .5);
+ fmdPlastic->DefineElement(1, 12.01, 6, .5);
+ fmdPlastic->SetFillColor(4);
+ fPlastic = new TGeoMedium("FMD Plastic", 1, fmdPlastic, pPlastic);
+
+ // Aluminium
+ Double_t pAl[] = { 0., mType, mMax, 10., .001, -1., .003, .003 };
+ TGeoMaterial* fmdAl = new TGeoMaterial("FMD Al", 26.981539, 13, 2.7);
+ fmdAl->SetFillColor(3);
+ fAl = new TGeoMedium("FMD Al", 1, fmdAl, pAl);
+
+ // Copper
+ Double_t pCopper[] = { 0., mType, mMax, 10., .01, 1., .003, .003 };
+ TGeoMaterial* fmdCopper = new TGeoMaterial("FMD Copper", 63.546, 29.,8.96);
+ fmdCopper->SetFillColor(3);
+ fCopper = new TGeoMedium("FMD Copper", 1, fmdCopper, pCopper);
+}
+
+#define DEGRAD TMath::Pi()/180.
+//____________________________________________________________________
+void
+Geometry::Detector2XYZ(UInt_t sector, UInt_t strip, TVector3& xyz)
+{
+ UInt_t mod = sector / 2;
+ if (!fMatricies) {
+ Warning("Detector2XYZ", "No matricies");
+ return;
+ }
+ TGeoMatrix* m = static_cast<TGeoMatrix*>(fMatricies->At(mod));
+ if (!m) {
+ Warning("Detector2XYZ", "No matrix found for module %d", mod);
+ return;
+ }
+ Debug(10, "Detector2XYZ", "Transforming (%2d,%3d)", sector, strip);
+ Double_t rmax = fB.Mod();
+ Double_t stripoff = fA.Mod();
+ Double_t dstrip = (rmax - stripoff) / fNStrips;
+ Double_t r = (strip + .5) * dstrip + stripoff; // fLowR
+ Double_t theta = ((sector % 2) - .5) * fTheta;
+ Double_t modThick = (fSiThickness
+ + fPrintboardThickness
+ + fCopperThickness
+ + fChipThickness
+ + fSpacer);
+ Debug(10,"Detector2XYZ", "Radius %7.3f, angle %7.3f (%f, %f)", r, theta,
+ fLowR, stripoff);
+ Double_t local[] = {
+ r * TMath::Cos(theta * DEGRAD),
+ r * TMath::Sin(theta * DEGRAD),
+ -modThick + fSiThickness / 2
+ };
+ Double_t master[3];
+ Debug(10, "Detector2XYZ", "Local (%7.3f,%7.3f,%7.3f)",
+ local[0], local[1], local[2]);
+ m->LocalToMaster(local, master);
+ Debug(10, "Detector2XYZ", "Master (%7.3f,%7.3f,%7.3f)",
+ master[0], master[1], master[2]);
+ xyz.SetXYZ(master[0], master[1], master[2]);
+}
+//____________________________________________________________________
+void
+Geometry::MakeDetailed(TGeoVolume* caveVolume)
+{
+ Info("MakeSimple", "Using a detailed geometry");
+ Double_t xv[6] = { fA.X(), fC.X(), fB.X(), fB.X(), fC.X(), fA.X() };
+ Double_t yv[6] = { fA.Y(), fC.Y(), fB.Y(), -fB.Y(), -fC.Y(), -fA.Y() };
+ Double_t rmax = fB.Mod();
+ Double_t stripoff = fA.Mod();
+ Double_t dstrip = (rmax - stripoff) / fNStrips;
+
+ // Double_t hybridThick = (fPrintboardThickness + fCopperThickness
+ // + fChipThickness);
+ // Double_t modThick = fSiThickness + fSpacer + hybridThick;
+ // Double_t fmdWidth = (modThick + fHoneycombThickness + fModuleSpacing
+ // + fLegLength + fSpacer);
+
+ // Top
+ // TGeoTube* fmdShape = new TGeoTube(fLowR-.1,fHighR+.1,fmdWidth/2+.1);
+ // TGeoPcon* fmdShape = new TGeoPcon(0, 360, 2);
+ // fmdShape->DefineSection(0, -fmdWidth / 2 - .1, fLowR-.1, fHighR+.1);
+ // fmdShape->DefineSection(1, fmdWidth / 2 + .1, fLowR-.1, fHighR+.1);
+ // TGeoVolume* fmdVolume = new TGeoVolume("FMD1", fmdShape, fVacuum);
+
+
+ // Sensor
+ TGeoXtru* sensorShape = new TGeoXtru(2);
+ sensorShape->DefinePolygon(6, xv, yv);
+ sensorShape->DefineSection(0, - fSiThickness/2);
+ sensorShape->DefineSection(1, fSiThickness/2);
+ TGeoVolume* sensorVolume= new TGeoVolume("FISE",sensorShape,fSi);
+ sensorVolume->SetLineColor(2);
+ // sensorVolume->VisibleDaughters(kFALSE);
+ // sensorVolume->SetVisibility(kTRUE);
+
+ // Virtual volume shape to divide
+ TGeoTubeSeg* activeShape = new TGeoTubeSeg(fLowR, rmax, fSiThickness/2,
+ - fTheta, fTheta);
+ TGeoVolume* activeVolume = new TGeoVolume("FIAC", activeShape,fSi);
+ activeVolume->SetLineColor(3);
+ sensorVolume->AddNodeOverlap(activeVolume, 0);
+ // activeVolume->SetTransparency(0x3f);
+ TGeoVolume* sectorVolume = activeVolume->Divide("FISC",2,2,-fTheta,
+ 0,0,"N");
+ TGeoVolume* stripVolume = sectorVolume->Divide("FIST",1,fNStrips,
+ stripoff,dstrip,0,"SX");
+ (void)stripVolume;
+
+ // Position
+ Double_t x, y, z;
+
+ // Make PCB volume
+ for (Int_t i = 0; i < 3; i++) yv[i] -= fBondingWidth;
+ for (Int_t i = 3; i < 6; i++) yv[i] += fBondingWidth;
+ Double_t off = (TMath::Tan(TMath::Pi() * fTheta / 180) * fBondingWidth);
+
+ // PCB layer
+ TGeoXtru* pcbShape = new TGeoXtru(2);
+ pcbShape->DefinePolygon(6, xv, yv);
+ pcbShape->DefineSection(0, - fPrintboardThickness/2);
+ pcbShape->DefineSection(1, fPrintboardThickness/2);
+ TGeoVolume* pcbVolume = new TGeoVolume("FPCB",pcbShape,fPCB);
+ pcbVolume->SetLineColor(4);
+
+ // Copper layer
+ TGeoXtru* cuShape = new TGeoXtru(2);
+ cuShape->DefinePolygon(6, xv, yv);
+ cuShape->DefineSection(0, - fCopperThickness/2);
+ cuShape->DefineSection(1, fCopperThickness/2);
+ TGeoVolume* cuVolume = new TGeoVolume("FCOP",cuShape,fCopper);
+ cuVolume->SetLineColor(4);
+
+ // Chip layer
+ TGeoXtru* chipShape = new TGeoXtru(2);
+ chipShape->DefinePolygon(6, xv, yv);
+ chipShape->DefineSection(0, - fChipThickness/2);
+ chipShape->DefineSection(1, fChipThickness/2);
+ TGeoVolume* chipVolume = new TGeoVolume("FCHI",chipShape,fChip);
+ chipVolume->SetLineColor(4);
+
+ // Short leg shape
+ TGeoTube* shortLegShape = new TGeoTube(0, fLegRadius, fLegLength / 2);
+ TGeoVolume* shortLegVolume = new TGeoVolume("FISL", shortLegShape, fPlastic);
+ shortLegVolume->SetLineColor(2);
+ // Long leg shape
+ TGeoTube* longLegShape = new TGeoTube(0, fLegRadius,
+ (fLegLength + fModuleSpacing) / 2);
+ TGeoVolume* longLegVolume = new TGeoVolume("FILL", longLegShape, fPlastic);
+ longLegVolume->SetLineColor(2);
+
+ // Make a front volume
+ TGeoVolume* frontVolume = new TGeoVolumeAssembly("FIFV");
+ z = fPrintboardThickness / 2;
+ frontVolume->AddNode(pcbVolume, 0, new TGeoTranslation(0, 0, z));
+ z += (fPrintboardThickness + fCopperThickness) / 2;
+ frontVolume->AddNode(cuVolume, 0, new TGeoTranslation(0, 0, z));
+ z += (fCopperThickness + fChipThickness) / 2;
+ frontVolume->AddNode(chipVolume, 0, new TGeoTranslation(0, 0, z));
+ z += (fChipThickness + fModuleSpacing
+ + fLegLength) / 2;
+ x = fA.X() + fLegOffset + fLegRadius;
+ y = 0;
+ frontVolume->AddNode(longLegVolume, 0, new TGeoTranslation(x, y, z));
+ x = fC.X();
+ y = fC.Y() - fLegOffset - fLegRadius - off;
+ frontVolume->AddNode(longLegVolume, 1, new TGeoTranslation(x,y,z));
+ y = -y;
+ frontVolume->AddNode(longLegVolume, 2, new TGeoTranslation(x,y,z));
+
+ // Make a back volume
+ TGeoVolume* backVolume = new TGeoVolumeAssembly("FIBV");
+ z = fPrintboardThickness / 2;
+ backVolume->AddNode(pcbVolume, 0, new TGeoTranslation(0, 0, z));
+ z += (fPrintboardThickness + fCopperThickness) / 2;
+ backVolume->AddNode(cuVolume, 0, new TGeoTranslation(0, 0, z));
+ z += (fCopperThickness + fChipThickness) / 2;
+ backVolume->AddNode(chipVolume, 0, new TGeoTranslation(0, 0, z));
+ z += (fChipThickness + fLegLength) / 2;
+ x = fA.X() + fLegOffset + fLegRadius;
+ y = 0;
+ backVolume->AddNode(shortLegVolume, 0, new TGeoTranslation(x, y, z));
+ x = fC.X();
+ y = fC.Y() - fLegOffset - fLegRadius - off;
+ backVolume->AddNode(shortLegVolume, 1, new TGeoTranslation(x,y,z));
+ y = -y;
+ backVolume->AddNode(shortLegVolume, 2, new TGeoTranslation(x,y,z));
+
+#if 1
+ TGeoVolume* ringTopVolume = new TGeoVolumeAssembly("FITV");
+ TGeoVolume* ringBotVolume = new TGeoVolumeAssembly("FIBV");
+ Double_t dz = 0;
+#else
+ Double_t dz = (fSiThickness + fSpacer + fModuleSpacing +
+ fPrintboardThickness + fCopperThickness +
+ fChipThickness + fLegLength);
+ TGeoTubeSeg* topRingShape = new TGeoTubeSeg(fA.X(), rmax, dz/2, 0, 180);
+ TGeoTubeSeg* botRingShape = new TGeoTubeSeg(fA.X(), rmax, dz/2, 180, 360);
+ TGeoVolume* ringTopVolume = new TGeoVolume("FITV", topRingShape, fAir);
+ TGeoVolume* ringBotVolume = new TGeoVolume("FIBV", botRingShape, fAir);
+#endif
+
+ TGeoVolume* half = ringTopVolume;
+ // Place in mother
+ for (Int_t i = 0; i < 10; i++) {
+ if (i > 4) half = ringBotVolume;
+ Bool_t front = ((i % 2) == 0);
+ z = -dz/2 + fSiThickness / 2 + (i % 2) * fModuleSpacing;
+ TGeoMatrix* m1 = new TGeoCombiTrans(0,0,z,0);
+ m1->RotateZ((2 * i + 1) * fTheta);
+ half->AddNode(sensorVolume, i, m1);
+ TGeoVolume* vol = (front ? frontVolume : backVolume);
+ z += fSpacer + fSiThickness / 2;
+ TGeoMatrix* m2 = new TGeoCombiTrans(0,0,z,0);
+ m2->RotateZ((2 * i + 1) * fTheta);
+ half->AddNode(vol, i, m2);
+ }
+
+ TGeoVolume* fmdTopVolume = new TGeoVolumeAssembly("FMT1");
+ TGeoVolume* fmdBotVolume = new TGeoVolumeAssembly("FMB1");
+ fmdTopVolume->AddNode(ringTopVolume, 0, new TGeoTranslation(0,0,dz/2));
+ fmdBotVolume->AddNode(ringBotVolume, 0, new TGeoTranslation(0,0,dz/2));
+
+ // Top of Honeycomb
+ TGeoTubeSeg* hcShape = new TGeoTubeSeg(fLowR, fHighR,
+ fHoneycombThickness / 2,
+ 0, 180);
+ TGeoVolume* hcVolume = new TGeoVolume("F1II", hcShape, fAl);
+ hcVolume->VisibleDaughters(kFALSE);
+ hcVolume->SetLineColor(5);
+ // Air in top of honeycomb
+ TGeoTubeSeg* ihcShape = new TGeoTubeSeg(fLowR + fAlThickness,
+ fHighR - fAlThickness,
+ (fHoneycombThickness
+ - fAlThickness) / 2,
+ 0, 180);
+ TGeoVolume* ihcVolume = new TGeoVolume("F1IK", ihcShape, fKapton);
+ hcVolume->AddNode(ihcVolume, 0);
+
+ // Add honey comb to mothers.
+ z = (fSiThickness + fSpacer + fPrintboardThickness + fCopperThickness
+ + fChipThickness + fModuleSpacing + fLegLength + fHoneycombThickness/2);
+ fmdTopVolume->AddNode(hcVolume, 0, new TGeoTranslation(0,0,z));
+ TGeoMatrix* r = new TGeoCombiTrans(0,0,z, 0); r->RotateZ(180);
+ fmdBotVolume->AddNode(hcVolume, 1, r);
+
+ z = 0;
+ caveVolume->AddNode(fmdTopVolume, 1, new TGeoTranslation(0,0,z));
+ caveVolume->AddNode(fmdBotVolume, 1, new TGeoTranslation(0,0,z));
+}
+
+//____________________________________________________________________
+void
+Geometry::Register()
+{
+ // Framework::Task::Register(option);
+ TGeoShape* caveShape = new TGeoBBox(fHighR+10, fHighR+10, fHighR+10);
+ TGeoVolume* caveVolume = new TGeoVolume("Cave", caveShape, fVacuum);
+ gGeoManager->SetTopVolume(caveVolume);
+
+ if (!gGeoManager->GetVolume("FM1T")) MakeDetailed(caveVolume);
+ gGeoManager->CloseGeometry();
+}
+
+//____________________________________________________________________
+void
+Geometry::Align()
+{
+ if (!gGeoManager) return;
+ if (!gGeoManager->GetTopVolume()) return;
+ if (!gGeoManager->IsClosed()) return;
+
+ if (!fMatricies) fMatricies = new TObjArray;
+ TGeoIterator next(gGeoManager->GetTopVolume());
+ TGeoNode* node = 0;
+
+ while ((node = static_cast<TGeoNode*>(next()))) {
+ // node->Print();
+ if (node->GetName()[0] == 'F' && node->GetName()[2] == 'S' &&
+ node->GetName()[3] == 'E') {
+ // We go an FMD module
+ TString path("/");
+ path.Append(gGeoManager->GetNode(0)->GetName());
+ Int_t nlevel = next.GetLevel();
+ Debug(10, "Exec", "Level is %d", nlevel);
+ for (int lvl = 0; lvl <= nlevel; lvl++) {
+ TGeoNode* p = next.GetNode(lvl);
+ if (!p) continue;
+ Debug(10, "Exec", "Adding '%s' to path '%s'",
+ p->GetName(), path.Data());
+ if (!path.IsNull()) path.Append("/");
+ path.Append(p->GetName());
+ }
+ TGeoPhysicalNode* pnode = gGeoManager->MakePhysicalNode(path.Data());
+ TGeoHMatrix* matrix = new TGeoHMatrix(*node->GetMatrix());
+ TGeoRotation pertub;
+ Double_t angles[] = { gRandom->Uniform(0, 3),
+ gRandom->Uniform(0, 3),
+ gRandom->Uniform(0, 3) };
+ pertub.RotateX(angles[0]);
+ pertub.RotateY(angles[1]);
+ pertub.RotateZ(angles[1]);
+ *matrix *= pertub;
+ Debug(5, "Exec", "Aliging %s (%f,%f,%f)",
+ pnode->GetName(), angles[0], angles[1], angles[2]);
+ pnode->Align(matrix);
+ }
+ }
+}
+
+//____________________________________________________________________
+//
+// EOF
+//
+
+
#include <TArrayI.h>
#include <AliRun.h>
#include <TNtuple.h>
-
+#include <AliDetector.h>
+#include <TFile.h>
+#include <TGeoManager.h>
+#include <TGeoVolume.h>
+#include <TGeoMedium.h>
+#include <TParticle.h>
+#include <THStack.h>
/** @class Media
@brief Media description
@ingroup FMD_script
*/
struct Media : public TNamed
{
- TArrayI* fMeds;
+ TArrayI fMeds;
TNtuple* fCount;
+ TH1D* fHist;
Media(const char* name, const char* title)
: TNamed(name, title), fMeds(0)
{
- fCount = new TNtuple(GetName(), "E/I:DET:RNG:SEC:STR:X/F:Y:Z:EDEP:PDG/I");
+ fHist = new TH1D(Form("h%s",name), title, 120, -6, 6);
+ // fHist->SetFillStyle(3001);
+ fCount = new TNtuple(GetName(), GetTitle(),
+ "E:DET:RNG:SEC:STR:X:Y:Z:EDEP:PDG");
}
- Media(const char* name)
- : TNamed(name, "Media information"), fMeds(0), fCount(0)
+ Media(const char* name, TArrayI* a)
+ : TNamed(name, "Media information"), fMeds(a->fN, a->fArray)
{
- AliDetector* det = gAlice->GetDetector(name);
- if (!det) {
- Warning("Media", "Detector %s not found in gAlice", name);
- return;
+ fHist = new TH1D(Form("h%s",name), GetTitle(), 120, -6, 6);
+ fHist->SetFillStyle(3001);
+ fHist->SetFillColor(fMeds.fN > 0 ? fMeds[0] : 0);
+ fCount = new TNtuple(GetName(), GetTitle(),
+ "E:DET:RNG:SEC:STR:X:Y:Z:EDEP:PDG");
+ std::cout << "Media list " << name << ":" << std::endl;
+ for (Int_t i = 0; i < fMeds.fN; i++) {
+ std::cout << " " << fMeds[i] << std::flush;
+ if (fMeds[i] == 0 && i > 0) break;
}
- fMeds = det->GetIdtmed();
- if (!fMeds) {
- Warning("Media", "No mediums for detector %s found", name);
- return;
+ std::cout << std::endl;
+ }
+ Bool_t HasMedia(Int_t id)
+ {
+ if (fMeds.fN <= 0) return kFALSE;
+ for (Int_t i = 0; i < fMeds.fN; i++) {
+ if (fMeds[i] == id) return kTRUE;
+ if (fMeds[i] == 0 && i > 0) break;
}
- fCount = new TNtuple(GetName(), "E/I:DET:RNG:SEC:STR:X/F:Y:Z:EDEP:PDG/I");
+ return kFALSE;
+ }
+ void Fill(Int_t ev, AliFMDHit* hit)
+ {
+ Float_t x[10];
+ x[0] = ev;
+ x[1] = hit->Detector();
+ x[2] = Int_t(hit->Ring());
+ x[3] = hit->Sector();
+ x[4] = hit->Strip();
+ x[5] = hit->X();
+ x[6] = hit->Y();
+ x[7] = hit->Z();
+ x[8] = hit->Edep();
+ x[9] = hit->Pdg();
+ // return;
+ fCount->Fill(x);
+ // r = TMath::Sqrt(hit->X() * hit->X() + hit->Y() * hit->Y());
+ // if(r!=0) Float_t wgt=1./(2. * 3.1415*r);
+ Double_t r = TMath::Sqrt(hit->X()*hit->X()+hit->Y()*hit->Y());
+ Double_t t = TMath::ATan2(r, hit->Z());
+ Double_t e = -TMath::Log(TMath::Tan(t/2));
+ fHist->Fill(e);
}
};
/** @class GetMedia
@brief Get media where a particle is produced
@code
+ aliroot
+ Root> gAlice->Init("$ALICE_ROOT/FMD/Config.C");
+ Root> .x $ALICE_ROOT/FMD/scriptsWriteMedArrays.C
+ Root> gAlice->Run();
+ Root> .q
+ aliroot
Root> .L Compile.C
Root> Compile("GetMedia.C")
Root> GetMedia c
TObjArray fMedia;
Media* fOther;
Media* fAll;
+ Media* fPrim;
Int_t fEv;
+ THStack* fStack;
TFile* fOutput;
public:
//__________________________________________________________________
GetMedia(const char* modlist="FMD:ITS:BODY:ABSO:START:PIPE",
const char* output="media.root")
+ : fModList(modlist)
{
+ AddLoad(kGeometry);
+ AddLoad(kKinematics);
+
fOutput = TFile::Open(output, "RECREATE");
- fOther = new Media("other", "Unknown media"),
- fAll = new Media("All", "All media")
+ fStack = new THStack("summed", "Summed hits");
+ fOther = AddMedia("other", "Unknown media");
+ fAll = AddMedia("All", "All media");
+ fPrim = AddMedia("Primary", "Primary particles");
fEv = 0;
- AddLoad(kKinematics);
+
+
+ fAll->fHist->SetFillColor(0);
+ fAll->fHist->SetFillStyle(3004);
+ fPrim->fHist->SetFillColor(1);
+ fPrim->fHist->SetFillStyle(3001);
+ fStack->Add(fPrim->fHist);
}
//__________________________________________________________________
- Media* FindMedia(Int_t med)
+ Media* AddMedia(const char* name, const char* title, TArrayI* a=0)
{
- TIter next(&fMedia);
Media* media = 0;
- while ((media == static_cast<Media*>(next()))) {
- if (!media->fMeds) continue;
- TArrayI& ids = *(media->fMeds);
- for (Int_t i = 0; i < ids.fN; i++)
- if (med == ids[i]) return media;
+ if (!a) media = new Media(name, title);
+ else media = new Media(name, a);
+ if (a) {
+ fMedia.Add(media);
+ fStack->Add(media->fHist);
+ }
+ return media;
+ }
+ //__________________________________________________________________
+ Media* FindMedia(Int_t med)
+ {
+ for (Int_t i = 0; i < fMedia.GetEntries(); i++) {
+ Media* media = static_cast<Media*>(fMedia.At(i));
+ if (!media) continue;
+ if (media->HasMedia(med)) return media;
}
return 0;
}
//__________________________________________________________________
Bool_t Init()
{
+ Bool_t ret = AliFMDInputHits::Init();
if (!gGeoManager) {
- Error("GetMedia", "No geometry defined - make sure you have that");
+ Error("Init", "No geometry defined - make sure you have that");
return kFALSE;
}
- if (!gAlice) {
- Error("GetMedia", "gAlice not defined");
+ if (!fRun) {
+ Error("Init", "gAlice not defined");
return kFALSE;
}
Int_t now = 0;
Int_t colon;
+ TFile* file = TFile::Open("medid.root", "READ");
+ if (!file) {
+ Error("Init", "couldn't open medid.root");
+ return kFALSE;
+ }
+ fOutput->cd();
while ((colon = fModList.Index(":", now)) >= 0) {
- fMedia.Add(new Media(fModList(now, colon-now)));
+ TString d(fModList(now, colon-now));
now = colon+1;
+ TArrayI* a = 0;
+ file->GetObject(d.Data(), a);
+ if (!a) {
+ Warning("Init", "No medium array for %s", d.Data());
+ continue;
+ }
+ AddMedia(d.Data(),0, a);
}
- if (now < fModList.Length())
- fMedia.Add(new Media(now, fModList.Length()-now));
+ if (now < fModList.Length()) {
+ colon = fModList.Length();
+ TString d(fModList(now, colon-now));
+ TArrayI* a = 0;
+ file->GetObject(d.Data(), a);
+ if (!a)
+ Warning("Init", "No medium array for %s", d.Data());
+ else
+ AddMedia(d.Data(),0, a);
+ }
+ file->Close();
if (fMedia.GetEntries() <= 0) return kFALSE;
- return AliFMDInputHits::Init();
+ return ret;
}
//__________________________________________________________________
/** Begining of event
<< track << ")" << std::endl;
return kFALSE;
}
+ fAll->Fill(fEv, hit);
+ if (track->IsPrimary()) {
+ fPrim->Fill(fEv, hit);
+ return kTRUE;
+ }
+
// Get production vertex
Double_t vx = track->Vx();
Double_t vy = track->Vy();
Double_t vz = track->Vz();
-
- fAll->fCounts->Fill(fEv, hit->Detector(),Int_t(hit->Ring()),
- hit->Sector(), hit->Strip(),
- hit->X(), hit->Y(), hit->Z(), hit->Edep(),
- hit->Pdg());
// Get node
TGeoNode* prodNode = gGeoManager->FindNode(vx,vy,vz);
- if (!prodNode) return kTRUE;
-
+ if (!prodNode) {
+ Warning("ProcessHit", "didn't find a node for production vertex");
+ return kTRUE;
+ }
+ //Info("ProcessHit", "Production vertex in node %s", prodNode->GetName());
+
// Get volume
TGeoVolume* prodVol = prodNode->GetVolume();
- if (!prodVol) return kTRUE;
+ if (!prodVol) {
+ Warning("ProcessHit", "didn't find a volume for production vertex");
+ return kTRUE;
+ }
+ //Info("ProcessHit", "Production vertex in volume %s",prodVol->GetName());
// Med medium
TGeoMedium* prodMed = prodVol->GetMedium();
- if (!prodMed) return kTRUE;
+ if (!prodMed) {
+ Warning("ProcessHit", "didn't find a medium for production vertex");
+ return kTRUE;
+ }
+ // Info("ProcessHit", "Production vertex in medium %s %d",
+ // prodMed->GetName(), prodMed->GetId());
- Media* media = FindMedia(prodMed->GetUniqueID());
- if (media) media = fOther;
- // r = TMath::Sqrt(hit->X() * hit->X() + hit->Y() * hit->Y());
- // if(r!=0) Float_t wgt=1./(2. * 3.1415*r);
- media->fCounts->Fill(fEv, hit->Detector(),Int_t(hit->Ring()),
- hit->Sector(), hit->Strip(),
- hit->X(), hit->Y(), hit->Z(), hit->Edep(),
- hit->Pdg());
+ Media* media = FindMedia(prodMed->GetId());
+ if (!media) {
+ Warning("ProcessHit", "Media not found %s (%d)",
+ prodMed->GetName(), prodMed->GetId());
+ media = fOther;
+ }
+ // Info("ProcessHit", "Adding to %s", media->GetName());
+ media->Fill(fEv, hit);
return kTRUE;
}
//__________________________________________________________________
Bool_t Finish()
{
+ fOutput->cd();
+ fStack->Write();
+ fStack->Draw();
fOutput->Write();
fOutput->Close();
+ fOutput = 0;
+ fMedia.Delete();
return kTRUE;
}
};
--- /dev/null
+//____________________________________________________________________
+//
+//
+// $Id$
+//
+// Script I used for rapid prototyping of the FMD3 geometry - in
+// particular the support cone
+//
+/** @defgroup node_geom Simple geometry
+ @ingroup FMD_script
+*/
+#include <TGeometry.h>
+#include <TNode.h>
+#include <TXTRU.h>
+#include <TTUBE.h>
+#include <TTUBS.h>
+#include <TPCON.h>
+#include <TBRIK.h>
+#include <TCanvas.h>
+#include <vector>
+#include <algorithm>
+#include <cmath>
+#include <iostream>
+
+//____________________________________________________________________
+/** @brief A 2D point
+ @ingroup node_geom
+ */
+struct point_t
+{
+ point_t(double x=0, double y=0) : first(x), second(y) {}
+ double first;
+ double second;
+};
+
+//____________________________________________________________________
+/** @brief Shape of a ring
+ @ingroup node_geom
+ */
+struct Ring
+{
+ // typedef std::pair<double,double> point_t;
+ typedef std::vector<point_t> points_t;
+ /** Constructor
+ @param rL Lower radius
+ @param rH Higer radius
+ @param theta Opening angle
+ @param waferR Wafer radius
+ @param siThick Silicon thickness
+ @param staggering Staggering of modules */
+ Ring(double rL, double rH, double theta, double waferR,
+ double siThick, double staggering)
+ : fStaggering(staggering),
+ fInnerRadius(rL),
+ fOuterRadius(rH),
+ fAngle(theta),
+ fRadius(waferR),
+ fThickness(siThick),
+ fVerticies(6)
+ {
+ double tan_theta = tan(fAngle * TMath::Pi() / 180.);
+ double tan_theta2 = pow(tan_theta,2);
+ double r2 = pow(fRadius,2);
+ double ir2 = pow(fInnerRadius,2);
+ double or2 = pow(fOuterRadius,2);
+ double y_A = tan_theta * fInnerRadius;
+ double x_D = fInnerRadius + sqrt(r2 - tan_theta2 * ir2);
+ double x_D2 = pow(x_D,2);
+ double y_B = sqrt(r2 - or2 + 2 * fOuterRadius * x_D - x_D2);
+ double x_C = ((x_D + sqrt(-tan_theta2 * x_D2
+ + r2 * (1 + tan_theta2)))
+ / (1 + tan_theta2));
+ double y_C = tan_theta * x_C;
+
+ fVerticies[0] = point_t(fInnerRadius, y_A);
+ fVerticies[1] = point_t(x_C, y_C);
+ fVerticies[2] = point_t(fOuterRadius, y_B);
+ fVerticies[3] = point_t(fOuterRadius, -y_B);
+ fVerticies[4] = point_t(x_C, -y_C);
+ fVerticies[5] = point_t(fInnerRadius, -y_A);
+ }
+ /** Destructor */
+ virtual ~Ring()
+ {
+ fVerticies.clear();
+ }
+ /** Create a shape
+ @return pointer to new shape */
+ TShape* CreateShape()
+ {
+ std::cout << "Creating Module shape" << std::flush;
+ TXTRU* moduleShape = new TXTRU("Module","Module", "", 6, 2);
+ for (Int_t i = 0; i < 6; i++) {
+ std::cout << "." << std::flush;
+ point_t& p = fVerticies[i];
+ moduleShape->DefineVertex(i, p.first, p.second);
+ }
+ moduleShape->DefineSection(0, -fThickness/2, 1, 0, 0);
+ moduleShape->DefineSection(1, fThickness/2, 1, 0, 0);
+ std::cout << std::endl;
+ return (TShape*)moduleShape;
+ }
+ /** Create a node that represents a ring.
+ @return Node */
+ TNode* CreateRing(const char* name, double z)
+ {
+ std::cout << "Creating Ring node for " << name << std::flush;
+ double bredth = fStaggering + fThickness;
+ TShape* ringShape = new TTUBE(Form("%sShape", name), "Ring Shape",
+ "", fInnerRadius,
+ fOuterRadius,bredth/2);
+ TNode* ringNode = new TNode(Form("%sNode", name), "Ring Node",
+ ringShape, 0, 0, z+bredth/2, 0);
+ TShape* moduleShape = CreateShape();
+ Int_t n = Int_t(360 / 2 / fAngle);
+ for (Int_t i = 0; i < n; i++) {
+ std::cout << "." << std::flush;
+ ringNode->cd();
+ Double_t theta = 2 * fAngle * i;
+ Double_t z = -(bredth+fThickness)/2+(i%2?0:fStaggering);
+ 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(2);
+ moduleNode->SetLineColor(2);
+ moduleNode->SetLineWidth(2);
+ }
+ std::cout << std::endl;
+ ringNode->SetVisibility(0);
+ return ringNode;
+ }
+ double fStaggering;
+ /** Inner radius */
+ double fInnerRadius;
+ /** Outer radius */
+ double fOuterRadius;
+ /** Opening angle (in degrees) */
+ double fAngle;
+ /** Radius (in centimeters) */
+ double fRadius;
+ /** Thickness */
+ double fThickness;
+ /** List of verticies */
+ points_t fVerticies;
+};
+
+//____________________________________________________________________
+/** @brief Shape of a detector
+ @ingroup node_geom
+ */
+struct Detector
+{
+ /** Constructor
+ @param id
+ @param inner
+ @param outer */
+ Detector(Ring* inner, double iZ, Ring* outer=0, double oZ=0)
+ : fInner(inner), fInnerZ(iZ), fOuter(outer), fOuterZ(oZ)
+ {}
+ /** Destructor */
+ virtual ~Detector() {}
+ /** Create rings */
+ virtual void CreateRings()
+ {
+ if (fInner) fInner->CreateRing("inner", fInnerZ);
+ if (fOuter) fOuter->CreateRing("outer", fOuterZ);
+ }
+ /** Create a node that represents the support */
+ virtual void CreateSupport(double) { }
+ /** Pointer to inner ring */
+ Ring* fInner;
+ /** Position in z of inner ring */
+ double fInnerZ;
+ /** Pointer to outer ring */
+ Ring* fOuter;
+ /** Position in z of inner ring */
+ double fOuterZ;
+};
+
+//____________________________________________________________________
+/** @brief FMD3 simple node geometry
+ @ingroup node_geom
+ */
+struct FMD3 : public Detector
+{
+ /** Constructor
+ @param inner Inner ring representation
+ @param outer Outer ring representation */
+ FMD3(Ring* inner, Ring* outer)
+ : Detector(inner, -62.8,outer, -75.2)
+ {
+ fNoseRl = 5.5;
+ fNoseRh = 6.7;
+ fNoseDz = 2.8 / 2;
+ fNoseZ = -46;
+ fConeL = 30.9;
+ fBackRl = 61 / 2;
+ fBackRh = 66.8 /2;
+ fBackDz = 1.4 / 2;
+ fBeamDz = .5 / 2;
+ fBeamW = 6;
+ fFlangeR = 49.25;
+ }
+ virtual ~FMD3() {}
+ void CreateRings()
+ {
+ double zdist = fConeL - 2 * fBackDz - 2 * fNoseDz;
+ double tdist = fBackRh - fNoseRh;
+ double alpha = tdist / zdist;
+ double x, rl, rh, z;
+ z = fNoseZ - fConeL / 2;
+ TPCON* fmd3Shape = new TPCON("fmd3Shape", "FMD 3 Shape", "", 0, 360, 7);
+ x = fNoseZ;
+ rl = fNoseRl;
+ rh = fNoseRh;
+ fmd3Shape->DefineSection(0, x - z, rl, rh);
+ x = fNoseZ-2*fNoseDz;
+ fmd3Shape->DefineSection(1, x - z, rl, rh);
+ x = fInnerZ - fInner->fStaggering - fInner->fThickness;
+ rl = fInner->fInnerRadius;
+ rh = fNoseRh + alpha * TMath::Abs(x-fNoseZ + 2 * fNoseDz);
+ fmd3Shape->DefineSection(2, x - z, rl, rh);
+ x = fOuterZ;
+ rl = fOuter->fInnerRadius;
+ rh = fBackRh;
+ fmd3Shape->DefineSection(3, x - z, rl, rh);
+ x = fNoseZ - zdist - 2 * fNoseDz;
+ rl = fOuter->fInnerRadius;
+ rh = fBackRh;
+ fmd3Shape->DefineSection(4, x - z, rl, rh);
+ x = fNoseZ - zdist - 2 * fNoseDz;
+ rl = fOuter->fInnerRadius;
+ rh = fFlangeR;
+ fmd3Shape->DefineSection(5, x - z, rl, rh);
+ x = fNoseZ - fConeL;
+ rl = fOuter->fInnerRadius;
+ rh = fFlangeR;
+ fmd3Shape->DefineSection(6, x - z, rl, rh);
+
+ TNode* fmd3Node = new TNode("fmd3Node", "FMD3 Node", fmd3Shape,0,0,z,0);
+ fmd3Node->SetLineColor(11);
+ fmd3Node->SetFillColor(11);
+ fmd3Node->SetVisibility(1);
+ fmd3Node->cd();
+ if (fInner) fInner->CreateRing("inner", fInnerZ-z);
+ fmd3Node->cd();
+ if (fOuter) fOuter->CreateRing("outer", fOuterZ-z);
+ fmd3Node->cd();
+ CreateSupport(fNoseZ - z);
+ }
+
+ /** Create support volumes */
+ void CreateSupport(double noseZ)
+ {
+ TShape* noseShape = new TTUBE("noseShape", "Nose Shape", "",
+ fNoseRl, fNoseRh, fNoseDz);
+ TNode* noseNode = new TNode("noseNode", "noseNode", noseShape,
+ 0, 0, noseZ - fNoseDz, 0);
+ noseNode->SetLineColor(0);
+ double zdist = fConeL - 2 * fBackDz - 2 * fNoseDz;
+ double tdist = fBackRh - fNoseRh;
+ double beamL = TMath::Sqrt(zdist * zdist + tdist * tdist);
+ double theta = -TMath::ATan2(tdist, zdist);
+ TShape* backShape = new TTUBE("backShape", "Back Shape", "",
+ fBackRl, fBackRh, fBackDz);
+ TNode* backNode = new TNode("backNode", "backNode", backShape,
+ 0, 0, noseZ-2*fNoseDz-zdist-fBackDz, 0);
+ backNode->SetLineColor(0);
+ TShape* beamShape = new TBRIK("beamShape", "beamShape", "",
+ fBeamDz, fBeamW / 2 , beamL / 2);
+ Int_t n = 8;
+ Double_t r = fNoseRl + 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*fNoseDz-zdist/2, beamRotation);
+ beamNode->SetLineColor(0);
+ }
+ Double_t flangel = (fFlangeR - fBackRh) / 2;
+ TShape* flangeShape = new TBRIK("flangeShape", "FlangeShape", "",
+ flangel, fBeamW / 2, fBackDz);
+ n = 4;
+ r = fBackRh + 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*fNoseDz-zdist-fBackDz,
+ flangeRotation);
+ flangeNode->SetLineColor(0);
+
+ }
+ }
+ /** Nose inner radius */
+ double fNoseRl;
+ /** Nose outer radius */
+ double fNoseRh;
+ /** Nose depth */
+ double fNoseDz;
+ /** Nose start position */
+ double fNoseZ;
+ /** Length of whole support structure */
+ double fConeL;
+ /** Inner radius of back ring */
+ double fBackRl;
+ /** Outer radius of back ring */
+ double fBackRh;
+ /** Thickness of back ring */
+ double fBackDz;
+ /** Thickness of beams */
+ double fBeamDz;
+ /** Width of beams */
+ double fBeamW;
+ /** Ending radius of flanges */
+ double fFlangeR;
+};
+
+//____________________________________________________________________
+/** @brief Create a node geometry
+ @ingroup node_geom
+ @code
+ .x NodeGeometry.C++
+ @endcode
+ */
+void
+NodeGeometry()
+{
+ 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();
+
+ Ring inner( 4.3, 17.2, 18, 13.4 / 2, .03, 1);
+ Ring outer(15.6, 28.0, 9, 13.4 / 2, .03, 1);
+ FMD3 fmd3(&inner, &outer);
+ fmd3.CreateRings();
+
+ TCanvas* c = new TCanvas("c", "c", 800, 800);
+ c->SetFillColor(1);
+ geometry->Draw();
+ // c->x3d("ogl");
+}
+//____________________________________________________________________
+//
+// EOF
+//
+++ /dev/null
-//____________________________________________________________________
-//
-//
-// $Id$
-//
-// Script I used for rapid prototyping of the FMD3 geometry - in
-// particular the support cone
-//
-/** @defgroup simple_geom Simple geometry
- @ingroup FMD_script
-*/
-/** Calculate wafer parameters
- @param dl Lower diameter
- @param dh Higer diameter
- @param theta Opening angle
- @param r Wafer radius
- @param verbose Maybe be verbose
- @return Corners of shape
- @ingroup simple_geom
-*/
-//____________________________________________________________________
-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;
-}
-
-/** Create the sensor shape
- @param name
- @param rl
- @param rh
- @param th
- @param r
- @param dz
- @ingroup simple_geom
- @return */
-//____________________________________________________________________
-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;
-}
-
-/**
- @param name
- @param rl
- @param rh
- @param th
- @param siThick
- @param waferR
- @param staggering
- @param z
- @ingroup simple_geom
- @return
-*/
-//____________________________________________________________________
-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;
-}
-
-//____________________________________________________________________
-/**
- @param noseRl
- @param noseRh
- @param noseDz
- @param noseZ
- @param backRl
- @param backRh
- @param backDz
- @param coneL
- @param beamW
- @param beamDz
- @param flangeR
- @return
- @ingroup simple_geom
-*/
-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;
-}
-
-
-
-
-//____________________________________________________________________
-/**
- @ingroup simple_geom
- */
-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");
-}
-//____________________________________________________________________
-//
-// EOF
-//
--- /dev/null
+//____________________________________________________________________
+//
+// $Id$
+//
+/** @file TestIndex.C
+ @author Christian Holm Christensen <cholm@nbi.dk>
+ @date Thu Mar 30 01:20:02 2006
+ @brief Test AliFMDIndex and AliFMDObjIndex
+*/
+#ifndef __CINT__
+#include <AliFMDIndex.h>
+#include <iostream>
+#include <TFile.h>
+#endif
+/** @defgroup index_test Test of AliFMDIndex and AliFMDObjIndex
+ @ingroup FMD_script
+*/
+/** Write an AliFMDIndex object to output stream
+ @ingroup index_test
+ @param o Stream
+ @param i Index object
+ @return @a o */
+std::ostream&
+operator << (std::ostream& o, const AliFMDIndex& i)
+{
+ return o << i.Name();
+}
+
+/** Do the comparision, and print to standard out
+ @ingroup index_test
+ @param lhs Left hand side
+ @param rhs Right hand side
+ @return @f$ lhs < rhs @f$ */
+bool
+cmp(const AliFMDIndex& lhs, const AliFMDIndex& rhs)
+{
+ bool ret = lhs < rhs;
+ std::cout << " (" << lhs << " < " << rhs << "): " << ret << std::endl;
+ return ret;
+}
+
+/** Test if two index objects are equivilant (are the same).
+ Equivilance is defined as
+ @f[
+ lhs \equiv rhs: \not (lhs < rhs \vee rhs < lhs)
+ @f]
+ @ingroup index_test
+ @param lhs Left hand side
+ @param rhs Right hand side
+ @return @c true if @a lhs and @a rhs are equivilant */
+bool
+equiv(const AliFMDIndex& lhs, const AliFMDIndex& rhs)
+{
+ bool ret = !(cmp(lhs,rhs) || cmp(rhs,lhs));
+ std::cout << " (" << lhs << " <> " << rhs << "): " << ret << std::endl;
+ return ret;
+}
+
+/** Check that @f$ \not (x < x)@f$
+ @ingroup index_test
+ @param x Object to test
+ @return @c true if @a x is not less than itself */
+bool
+self(const AliFMDIndex& x)
+{
+ bool ret = !cmp(x,x);
+ std::cout << " !(" << x << " < " << x << "): " << ret << std::endl;
+ return ret;
+}
+
+/** Check if @f$ a \wedge b \Rightarrow c@f$
+ @ingroup index_test
+ @param a First condition
+ @param b Second condition
+ @param c Implication
+ @return @c true if the implication is valid */
+bool
+imply(bool a, bool b, bool c)
+{
+ bool ret = ((a && b) && c) || (!(a && b));
+ return ret;
+}
+
+/** Check if the comparison operator is transitive, that is
+ @f[
+ (x < y \wedge y < z) \Rightarrow x < z
+ @f]
+ @ingroup index_test
+ @param x First object
+ @param y Second object
+ @param z Third object
+ @return @c true if the implication is met. */
+bool
+trans(const AliFMDIndex& x, const AliFMDIndex& y, const AliFMDIndex& z)
+{
+ bool ret = imply(cmp(x,y), cmp(y,z), cmp(x,z));
+ std::cout << " (" << x << " < " << y << " && " << y << " < " << z
+ << ") => " << x << " < " << z << " "
+ << (ret ? "holds" : "violated") << std::endl;
+ return ret;
+
+}
+
+/** Check that the comparison operator preserves equivilance, that is
+ @f[
+ (x \equiv y \wedge y \equiv z) \Rightarrow (x \equiv z)
+ @f]
+ @ingroup index_test
+ @param x First object
+ @param y Second argument
+ @param z Third object
+ @return @c true if the implication holds. */
+bool
+equiv(const AliFMDIndex& x, const AliFMDIndex& y, const AliFMDIndex& z)
+{
+ bool ret = imply(equiv(x,y), equiv(y,z), equiv(x,z));
+ std::cout << " (" << x << " <> " << y << " && " << y << " <> " << z
+ << ") => " << x << " <> " << z << " "
+ << (ret ? "holds" : "violated") << std::endl;
+ return ret;
+}
+
+
+/** Check if the comparison operator is a @e strictly @e weak @e ordering
+ @ingroup index_test
+ */
+void
+TestIndex()
+{
+ AliFMDIndex i1(1, 'I', 5, 63);
+ AliFMDIndex i2(i1);
+ AliFMDIndex i3(i1);
+ AliFMDIndex i4(1, 'I', 5, 127);
+ AliFMDIndex i5(1, 'I', 15, 63);
+ AliFMDIndex i6(2, 'O', 15, 60);
+ std::cout << "Is !(x < x): " << std::endl;
+ self(i1);
+ std::cout << "Does (x < y && y < z) imply x < z: " << std::endl;
+ // true, true
+ trans(i1, i4, i5);
+ // true, false
+ trans(i1, i6, i5);
+ // false, true
+ trans(i4, i1, i5);
+ // false, false
+ trans(i6, i5, i1);
+ std::cout << "Does !(x < y || x > y) && !(y < z || z > y) imply "
+ << "!(x < z || x > z)" << std::endl;
+ // true, true
+ equiv(i1, i2, i3);
+ // true, false
+ equiv(i1, i2, i4);
+ // false, true
+ equiv(i4, i1, i2);
+ // false, false
+ equiv(i1, i4, i5);
+
+
+
+ TFile* file = TFile::Open("index.root", "RECREATE");
+ file->WriteObject(&i1,"i1");
+ file->Write();
+ file->Close();
+
+ AliFMDIndex* i7 = 0;
+ file = TFile::Open("index.root", "READ");
+ file->GetObject("i1", i7);
+ file->Close();
+ std::cout << *i7 << " == " << i1 << ": " << (*i7 == i1) << std::endl;
+
+}
+
+/** Check if the comparison operator is a @e strictly @e weak @e ordering
+ @ingroup index_test
+ */
+void
+TestObjIndex()
+{
+ AliFMDObjIndex i1(1, 'I', 5, 63);
+ AliFMDObjIndex i2(i1);
+ AliFMDObjIndex i3(i1);
+ AliFMDObjIndex i4(1, 'I', 5, 127);
+ AliFMDObjIndex i5(1, 'I', 15, 63);
+ AliFMDObjIndex i6(2, 'O', 15, 60);
+ std::cout << "Is !(x < x): " << std::endl;
+ self(i1);
+ std::cout << "Does (x < y && y < z) imply x < z: " << std::endl;
+ // true, true
+ trans(i1, i4, i5);
+ // true, false
+ trans(i1, i6, i5);
+ // false, true
+ trans(i4, i1, i5);
+ // false, false
+ trans(i6, i5, i1);
+ std::cout << "Does !(x < y || x > y) && !(y < z || z > y) imply "
+ << "!(x < z || x > z)" << std::endl;
+ // true, true
+ equiv(i1, i2, i3);
+ // true, false
+ equiv(i1, i2, i4);
+ // false, true
+ equiv(i4, i1, i2);
+ // false, false
+ equiv(i1, i4, i5);
+
+
+
+ TFile* file = TFile::Open("index.root", "RECREATE");
+ i1.Write("i1");
+ file->Write();
+ file->Close();
+
+ file = TFile::Open("index.root", "READ");
+ AliFMDObjIndex* i7 = (AliFMDObjIndex*)(file->Get("i1"));
+ file->Close();
+ std::cout << *i7 << " == " << i1 << ": " << (*i7 == i1) << std::endl;
+
+}
+
+/** Check that we can sort an array of index objects
+ @ingroup index_test
+ */
+void
+SortIndex()
+{
+ TList l;
+ for (size_t i = 0; i < 30; i++) {
+ UShort_t det = gRandom->Integer(3)+1;
+ Char_t ring = (gRandom->Uniform() > .5 ? 'O' : 'I');
+ UShort_t sec = gRandom->Integer(ring == 'I' ? 20 : 40);
+ UShort_t str = gRandom->Integer(ring == 'I' ? 512 : 256);
+ l.AddAt(new AliFMDObjIndex(det, ring, sec, str), i);
+ }
+ std::cout << "Before sort" << std::endl;
+ l.ls();
+ l.Sort();
+ std::cout << "After sort" << std::endl;
+ l.ls();
+}
+
+
+//____________________________________________________________________
+//
+// EOF
+//
--- /dev/null
+void
+WriteMedArrays()
+{
+ TFile* file = TFile::Open("medid.root", "RECREATE");
+ if (!file) {
+ Warning("WriteMedArrays", "failed to open medid.root");
+ return;
+ }
+ TObjArray* modules = gAlice->Modules();
+ if (!modules) {
+ Warning("WriteMedArrays", "failed to get modules");
+ return;
+ }
+ TIter next(modules);
+ AliModule* module = 0;
+ while ((module = static_cast<AliModule*>(next()))) {
+ Info("WriteMedArrays", "Getting medium id's for %s", module->GetName());
+ TArrayI* mediumIds = module->GetIdtmed();
+ if (!mediumIds) {
+ Warning("WriteMedArrays", "No medium id's for %s", module->GetName());
+ continue;
+ }
+ file->WriteObject(mediumIds,module->GetName());
+ }
+ file->Write();
+ file->Close();
+}
+
+
+