+++ /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$ */
-
-//
-// Realisation of AliLhcProcess for the fast simulation of the
-// Intra Beam Scattering process
-// in transverse and longitudinal direction.
-// Author: Andreas Morsch
-// andreas.morsch@cern.ch
-//
-
-#include <TCanvas.h>
-#include <TF1.h>
-#include <TGraph.h>
-#include <TH1F.h>
-#include <TMath.h>
-#include <TMultiGraph.h>
-
-#include "AliLhcProcessIBS.h"
-#include "AliLHC.h"
-#include "AliLhcIRegion.h"
-#include "AliLhcBeam.h"
-
-ClassImp(AliLhcProcessIBS)
-
- Double_t func(Double_t *x, Double_t *par);
-
-AliLhcProcessIBS::AliLhcProcessIBS(AliLHC* lhc, const char* name, const char* title)
- :AliLhcProcess(lhc,name,title),
- fCrossSection(0.),
- fIRegions(0),
- fTaux(0.),
- fTaue(0.),
- fTauxArray(0),
- fTaueArray(0)
-{
-// Constructor
- for (Int_t i = 0; i < 2; i++) {
- fBeam[i] = 0;
- fR[i] = 0.;
- fE[i] = 0.;
- }
-}
-
-AliLhcProcessIBS::AliLhcProcessIBS(const AliLhcProcessIBS& ibs):
- AliLhcProcess(ibs),
- fCrossSection(0.),
- fIRegions(0),
- fTaux(0.),
- fTaue(0.),
- fTauxArray(0),
- fTaueArray(0)
-{
-// Copy Constructor
- for (Int_t i = 0; i < 2; i++) {
- fBeam[i] = 0;
- fR[i] = 0.;
- fE[i] = 0.;
- }
-}
-
-
-AliLhcProcessIBS::~AliLhcProcessIBS()
-{
-// Destructor
-
-}
-
-void AliLhcProcessIBS::Init()
-{
- // Initialization
- const Float_t r0=1.535e-16;
-
- printf("\n Initializing Process %s", GetName());
- printf("\n ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^");
-
- fIRegions = fAccelerator->IRegions();
-
- for (Int_t i = 0; i < 2; i++) {
- fBeam[i] = fAccelerator->Beam(i);
- fR[i] = r0*fBeam[i]->Z()*fBeam[i]->Z()/fBeam[i]->A();
- fE[i] = 0.938*fBeam[i]->A();
- }
-}
-
-void AliLhcProcessIBS::Evolve(Float_t dt)
-{
-//
-// Evolve by one time step dt
- printf("\n Here process %s %f:", GetName(), dt);
- for (Int_t i=0; i<2; i++) {
- // Density
- Float_t sige = fBeam[i]->EnergySpread();
- Float_t avbeta = fAccelerator->AverageBeta();
- Float_t avd = fAccelerator->AverageDisp();
-
- Float_t gamma = fBeam[i]->Gamma();
- Float_t epsx = fBeam[i]->Emittance()*gamma;
- Float_t epsy = epsx;
- Float_t epse = fBeam[i]->LongEmittance();
- Float_t sigxb = TMath::Sqrt(epsx/gamma*avbeta);
- Float_t sigx = TMath::Sqrt(sigxb*sigxb+(avd*avd*sige*sige));
- Float_t ssigx = TMath::Sqrt(epsx/gamma/avbeta);
- Float_t sigy = sigx;
- Float_t ssigy = ssigx;
-
- Float_t asd = fBeam[i]->N()*fR[i]*fR[i]*fE[i]/
- (16.*TMath::Pi()*gamma*epsx*epsy*epse);
-
- // impact parameter
- Float_t d = sige*avd/sigx;
- //
- Float_t at = sige*TMath::Sqrt((1.-d)*(1.+d))/(gamma*ssigx);
- Float_t bt = sige*TMath::Sqrt((1.-d)*(1.+d))/(gamma*ssigy);
- Float_t ct = sige*TMath::Sqrt((1.-d)*(1.+d))*TMath::Sqrt(4.*sigy/fR[i]);
- //
- Double_t par[3];
- par[0] = at;
- par[1] = bt;
- par[2] = ct;
- TF1 *fct = new TF1("func",func, 0., 1., 3);
-
- Float_t f = (Float_t) 8.0*TMath::Pi()*
- fct->Integral(0., 1., par, 1.e-5);
- //
- const Double_t osq2=1./TMath::Sqrt(2.);
- fTaux = 1./(asd*f*(d-osq2*at)*(d+osq2*at));
- fTaue = 1./(asd*f*((1.-d)*(1.+d)));
- // printf("\n taux, taue %f %f", taux, taue);
- // Float_t tauy = -2./at*at/asd/f;
- fBeam[i]->IncreaseEmittance(dt/fTaux, dt/fTaue);
- }
-
-}
-
-void AliLhcProcessIBS::SetMonitor(Int_t n)
-{
- // Initialize Monitor
- if (fTauxArray) delete fTauxArray;
- if (fTaueArray) delete fTaueArray;
- fTauxArray = new Float_t[n];
- fTaueArray = new Float_t[n];
- fNmax = n;
-}
-
-void AliLhcProcessIBS::Record()
-{
- // Record monitor quantities
- fTauxArray[fAccelerator->Nt()] = fTaux/3600.;
- fTaueArray[fAccelerator->Nt()] = fTaue/3600.;
-}
-
-
-void AliLhcProcessIBS::DrawPlots()
-{
- // Draw monitor plots
- Float_t* t = fAccelerator->TimeA();
-
- TH1 *t1 = new TH1F("t1","Hor. IBS growth time",fNmax,0,t[fNmax]);
- t1->SetMinimum(0);
- t1->SetMaximum(fTauxArray[fNmax]*1.1);
- t1->SetStats(0);
- t1->GetXaxis()->SetTitle("t (h)");
- t1->GetYaxis()->SetTitle("tau_x (t)");
-
- TH1 *t2 = new TH1F("t2","Long. IBS growth time",fNmax,0,t[fNmax]);
- t2->SetMinimum(0);
- t2->SetMaximum(fTaueArray[fNmax]*1.1);
- t2->SetStats(0);
- t2->GetXaxis()->SetTitle("t (h)");
- t2->GetYaxis()->SetTitle("tau_l (t)");
-
- TGraph* grTaux = new TGraph(fNmax, fAccelerator->TimeA(), fTauxArray);
- grTaux->SetHistogram(t1);
-
- TGraph* grTaue = new TGraph(fNmax, fAccelerator->TimeA(), fTaueArray);
- grTaue->SetHistogram(t2);
- grTaue->SetLineStyle(2);
-
- TMultiGraph* mg = new TMultiGraph();
- mg->Add(grTaux);
- mg->Add(grTaue);
-
- TCanvas *c3 = new TCanvas("c3","IBS", 200, 10, 700, 500);
- c3->SetGrid();
- mg->Draw("AC");
- mg->GetXaxis()->SetTitle("t (h)");
- mg->GetYaxis()->SetTitle("IBS Growth Time (h)");
- mg->Draw("AC");
-}
-
-
-
-
-
-AliLhcProcessIBS& AliLhcProcessIBS::operator=(const AliLhcProcessIBS & /*rhs*/)
-{
-// Assignment operator
- return *this;
-}
-
-Double_t func(Double_t *x, Double_t *par)
-{
- Double_t a = par[0];
- Double_t b = par[1];
- Double_t cc = par[2];
-
- const Double_t kbc = 0.5772;
- Double_t xx = x[0];
- Double_t x2=xx*xx;
- Double_t x1=1.0-x2;
- Double_t p=1.0/TMath::Sqrt(x2+a*a*x1);
- Double_t q=1.0/TMath::Sqrt(x2+b*b*x1);
- return (1.0-3.0*x2)*p*q*(2.0*TMath::Log(0.5*cc*(p+q))-kbc);
-}
-
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff