Fix for the loophole in the magnets currents check: the L3Off/DipON was passing the...
[u/mrichter/AliRoot.git] / LHC / AliLhcIRegion.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 //
19 // Realisation of AliLhcMonitor simulating an LHC interaction region.
20 // The interaction region is described by the two LHC beams,
21 // by the beta* and the crossing angle. 
22 // As a monitor it records the luminosity, average luminosity and beta*
23 // time evolution.
24 // Author: Andreas Morsch
25 // andreas.morsch@cern.ch
26 //
27 #include "AliLhcIRegion.h"
28 #include "AliLhcBeam.h"
29 #include "AliLHC.h"
30
31 #include <TMath.h>
32 #include <TCanvas.h>
33 #include <TPad.h>
34 #include <TGraph.h>
35 #include <TMultiGraph.h>
36 #include <TH1.h>
37
38 ClassImp(AliLhcIRegion)
39
40 AliLhcIRegion::AliLhcIRegion(AliLHC* lhc, const char* name, const char* title)
41     :TNamed(name,title),
42      fAccelerator(lhc),
43      fBeam1(0),
44      fBeam2(0),
45      fLuminosity(0.),
46      fLuminosity0(0.),
47      fAverageLumi(0.),
48      fBetaStar(0.),
49      fBetaStar0(0.),
50      fCrossingAngle(0.),
51      fFrequency(0.),
52      fLumiArray(0),
53      fAverageLumiArray(0),
54      fBetaStarArray(0)
55 {
56 // Constructor
57 }
58
59 AliLhcIRegion::AliLhcIRegion(const AliLhcIRegion& region): 
60     TNamed(region), AliLhcMonitor(region),
61     fAccelerator(0),
62     fBeam1(0),
63     fBeam2(0),
64     fLuminosity(0.),
65     fLuminosity0(0.),
66     fAverageLumi(0.),
67     fBetaStar(0.),
68     fBetaStar0(0.),
69     fCrossingAngle(0.),
70     fFrequency(0.),
71     fLumiArray(0),
72     fAverageLumiArray(0),
73     fBetaStarArray(0)
74 {
75 // copy constructor
76 }
77
78 AliLhcIRegion::~AliLhcIRegion()
79 {
80 // Destructor
81
82   if (fLumiArray)        delete fLumiArray;
83   if (fAverageLumiArray) delete fAverageLumiArray;
84   if (fBetaStarArray)    delete fBetaStarArray;
85 }
86
87
88 AliLhcIRegion& AliLhcIRegion::operator=(const  AliLhcIRegion & /*rhs*/)
89 {
90 // Assignment operator
91     return *this;
92 }
93
94 void AliLhcIRegion::Init()
95 {
96   // Initialization
97   printf("\n Initializing Interaction Region %s", GetName());
98   printf("\n ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^");
99   // Initial Luminosity
100   fBeam1 = fAccelerator->Beam(0);
101   fBeam2 = fAccelerator->Beam(1);
102   fFrequency = 3.e10/(2.*TMath::Pi()*fAccelerator->Radius());
103
104   Luminosity();
105
106   fLuminosity0 = fLuminosity;  
107   fAverageLumi=0.;
108   fBetaStar0   = fBetaStar;
109   printf("\n IR Beta*                 :%10.3e cm        ", fBetaStar);
110   printf("\n IR Initial Luminosity    :%10.3e cm^-2s^-1 ", fLuminosity);
111 }
112
113 Float_t AliLhcIRegion::Luminosity()
114 {
115   Float_t sigma1 = TMath::Sqrt(fBeam1->Emittance()*fBetaStar);
116   Float_t sigma2 = TMath::Sqrt(fBeam2->Emittance()*fBetaStar);
117   fLuminosity = fFrequency * 
118     fBeam1->N()*fBeam2->N()/(2.*TMath::Pi()*(sigma1*sigma1+sigma2*sigma2));
119   return fLuminosity;
120 }
121
122 void AliLhcIRegion::Update()
123 {
124   Luminosity();
125 }
126
127 void AliLhcIRegion::SetMonitor(Int_t n)
128 {
129   // Initialize monitors
130
131   if (fLumiArray)        delete fLumiArray;
132   if (fAverageLumiArray) delete fAverageLumiArray;
133   if (fBetaStarArray)    delete fBetaStarArray;
134
135   fLumiArray        = new Float_t[n];
136   fAverageLumiArray = new Float_t[n];
137   fBetaStarArray    = new Float_t[n];
138
139   fAverageLumiArray[0] = 0.;
140   fNmax = n;
141 }
142
143 void AliLhcIRegion::Record()
144 {
145   // Record some time dependent quantities
146   //
147   Int_t n = fAccelerator->Nt();
148   // Luminosity
149   
150   fLumiArray[n] = fLuminosity;
151   
152   // Average Luminosity respecting set-up and filling time
153   if (fAccelerator->Time() > fAccelerator->SetUpTime())
154     fAverageLumi+=fLuminosity*fAccelerator->TimeStep();
155   
156   fAverageLumiArray[n] = fAverageLumi/
157     (Float_t(n+1)*fAccelerator->TimeStep()+fAccelerator->FillingTime())/
158     fLuminosity0;
159   
160   // Beta* 
161   fBetaStarArray[n] = fBetaStar;
162 }
163
164
165 void AliLhcIRegion::DrawPlots()
166 {
167   //
168   // Draw the monitor plots
169   //
170   Float_t* t = fAccelerator->TimeA();
171   //
172   char name1[20], name2[20], hname[20];
173   sprintf(name1,"c%s",GetName());
174   sprintf(name2,"b%s",GetName());
175   char title[30];
176   sprintf(title,"Luminosity Lifetime for %s",GetName());
177   
178   //
179   sprintf(hname,"%s%d",name1,0);
180   TH1 *g1 = new TH1F(hname,"Luminosity",fNmax,0,t[fNmax]);
181   g1->SetMinimum(0);
182   g1->SetMaximum(fLumiArray[0]*1.1);
183   g1->SetStats(0);
184   g1->GetXaxis()->SetTitle("t (h)");
185   g1->GetYaxis()->SetTitle("L(t) (cm**-2 s**-1)");
186   sprintf(hname,"%s%d",name1,1);
187   TH1 *g2 = new TH1F(hname,"Luminosity",fNmax,0,t[fNmax]);
188   g2->SetMinimum(0);
189   g2->SetMaximum(1.1);
190   g2->SetStats(0);
191   g2->GetXaxis()->SetTitle("t (h)");
192   g2->GetYaxis()->SetTitle("L(t)/L0");
193   sprintf(hname,"%s%d",name1,3);
194
195   TH1 *g3 = new TH1F(hname,"Average Luminosity",fNmax,0,t[fNmax]);
196   g3->SetMinimum(0);
197   g3->SetMaximum(1.1);
198   g3->SetStats(0);
199   
200   g3->GetXaxis()->SetTitle("t (h)");
201   g3->GetYaxis()->SetTitle("L(t)/L0");
202   sprintf(hname,"%s%d",name1,3);
203
204   TH1 *g4 = new TH1F(hname,"Beta*",fNmax,0,t[fNmax]);
205   g4->SetMinimum(0);
206   g4->SetMaximum(fBetaStarArray[0]*1.1);
207   g4->SetStats(0);
208   g4->GetXaxis()->SetTitle("t (h)");
209   g4->GetYaxis()->SetTitle("beta* (cm)");
210
211   TGraph* grLumi = new TGraph(fNmax, t, fLumiArray);
212   grLumi->SetHistogram(g1);
213
214   for (Int_t i=0; i<fNmax; i++) {
215     fLumiArray[i]=fLumiArray[i]/fLuminosity0;
216   }
217   TGraph* grLumiN = new TGraph(fNmax, t, fLumiArray);
218   grLumiN->SetHistogram(g2);
219
220   TGraph* grLumiA = new TGraph(fNmax, t, fAverageLumiArray);
221   grLumiA->SetHistogram(g3);
222   grLumiA->SetLineStyle(2);
223
224   TGraph* grLumiB = new TGraph(fNmax, t, fBetaStarArray);
225   grLumiB->SetHistogram(g4);
226
227   TMultiGraph* mg = new TMultiGraph();
228   mg->Add(grLumiN);
229   mg->Add(grLumiA);
230
231  
232
233   TCanvas *c1 = new TCanvas(name1,title, 200, 10, 700, 500);
234   
235   c1->SetGrid();
236   mg->Draw("AC");
237   mg->GetXaxis()->SetTitle("t (h)");
238   mg->GetYaxis()->SetTitle("L(t)/L0 and <L>/L0");
239   mg->Draw("AC");
240
241   TCanvas *c2 = new TCanvas(name2,title, 200, 10, 700, 500);
242   c2->SetGrid();
243   grLumiB->Draw("AC");
244
245 }
246
247
248
249