STEER/AliMagFC.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 //     Constant magnetic field class
  20 //     Used by AliRun class
  21 //     Author:
  22 //-------------------------------------------------------------------------
  23
  24 #include <stdlib.h>
  25
  26 #include "AliLog.h"
  27 #include "AliMagFC.h"
  28
  29 ClassImp(AliMagFC)
  30
  31 //________________________________________
  32 AliMagFC::AliMagFC()
  33     :AliMagF(),
  34     fCompensator(kFALSE),
  35     fBeamType(kBeamTypepp),
  36     fBeamEnergy(0),
  37     fQuadGradient(0),
  38     fDipoleField(0),
  39     fCCorrField(0),
  40     fACorr1Field(0),
  41     fACorr2Field(0)
  42 {
  43   //
  44   // Default constructor
  45   //
  46 }
  47
  48 //________________________________________
  49 AliMagFC::AliMagFC(const char *name, const char *title, Int_t integ,
  50                    Float_t factor, Float_t fmax)
  51     : AliMagF(name,title,integ,factor,fmax),
  52     fCompensator(kFALSE),
  53     fBeamType(kBeamTypepp),
  54     fBeamEnergy(7000.),
  55     fQuadGradient(0),
  56     fDipoleField(0),
  57     fCCorrField(0),
  58     fACorr1Field(0),
  59     fACorr2Field(0)
  60
  61 {
  62   //
  63   // Standard constructor
  64   //
  65   fType = kConst;
  66   fMap  = 1;
  67
  68 }
  69
  70 //________________________________________
  71 void AliMagFC::Field(Float_t *x, Float_t *b) const
  72 {
  73   //
  74   // Method to return the field in a point
  75   //
  76   b[0]=b[1]=b[2]=0;
  77   if(fMap==1) {
  78     if(TMath::Abs(x[2])<700 && x[0]*x[0]+(x[1]+30)*(x[1]+30) < 560*560) {
  79       b[2]=2;
  80     }
  81     else {
  82       if(-725 >= x[2] && x[2] >= -1225 ){
  83         Float_t dz = TMath::Abs(-975-x[2])*0.01;
  84         b[0] = - (1-0.1*dz*dz)*7;
  85         if(fFactor!=1) {
  86             b[0]*=fFactor;
  87             b[1]*=fFactor;
  88             b[2]*=fFactor;
  89         }
  90       }
  91       else {
  92           ZDCField(x, b);
  93       }
  94     }
  95
  96   }
  97   else {
  98       AliFatal(Form("Invalid field map for constant field %d",fMap));
  99   }
 100 }
 101
 102 //___________________________________________________
 103 void AliMagFC::ZDCField(Float_t *x, Float_t *b) const
 104 {
 105   // ---- This is the ZDC part
 106
 107   Float_t rad2 = x[0] * x[0] + x[1] * x[1];
 108   static Bool_t init = kFALSE;
 109
 110   if (! init) {
 111       init = kTRUE;
 112       //////////////////////////////////////////////////////////////////////
 113       // ---- Magnetic field values (according to beam type and energy) ----
 114       if(fBeamType==kBeamTypepp && fBeamEnergy == 5000.){
 115           // p-p @ 5+5 TeV
 116           fQuadGradient = 15.7145;
 117           fDipoleField  = 27.0558;
 118           // SIDE C
 119           fCCorrField   = 9.7017;
 120           // SIDE A
 121           fACorr1Field  = -13.2143;
 122           fACorr2Field  = -11.9909;
 123       } else if (fBeamType == kBeamTypepp && fBeamEnergy == 450.) {
 124           // p-p 0.45+0.45 TeV
 125           Float_t const kEnergyRatio = fBeamEnergy / 7000.;
 126
 127           fQuadGradient = 22.0002 * kEnergyRatio;
 128           fDipoleField  = 37.8781 * kEnergyRatio;
 129           // SIDE C
 130           fCCorrField   =  9.6908;
 131           // SIDE A
 132           fACorr1Field  = -13.2014;
 133           fACorr2Field  = -9.6908;
 134       } else if ((fBeamType == kBeamTypepp && fBeamEnergy == 7000.) ||
 135                  (fBeamType == kBeamTypeAA))
 136       {
 137           // Pb-Pb @ 2.7+2.7 TeV or p-p @ 7+7 TeV
 138           fQuadGradient = 22.0002;
 139           fDipoleField  = 37.8781;
 140           // SIDE C
 141           fCCorrField   = 9.6908;
 142           // SIDE A
 143           fACorr1Field  = -13.2014;
 144           fACorr2Field  = -9.6908;
 145       }
 146   }
 147
 148
 149   // SIDE C **************************************************
 150   if(x[2]<0.){
 151     if(x[2] < kCCorrBegin && x[2] > kCCorrEnd && rad2 < kCCorrSqRadius){
 152         if (fFactor != 0.) {
 153             b[0] = fCCorrField;
 154             b[1] = 0.;
 155             b[2] = 0.;
 156         }
 157     }
 158     else if(x[2] < kCQ1Begin && x[2] > kCQ1End && rad2 < kCQ1SqRadius){
 159         b[0] = fQuadGradient*x[1];
 160         b[1] = fQuadGradient*x[0];
 161         b[2] = 0.;
 162     }
 163     else if(x[2] < kCQ2Begin && x[2] > kCQ2End && rad2 < kCQ2SqRadius){
 164         b[0] = -fQuadGradient*x[1];
 165         b[1] = -fQuadGradient*x[0];
 166         b[2] = 0.;
 167     }
 168     else if(x[2] < kCQ3Begin && x[2] > kCQ3End && rad2 < kCQ3SqRadius){
 169         b[0] = -fQuadGradient*x[1];
 170         b[1] = -fQuadGradient*x[0];
 171         b[2] = 0.;
 172     }
 173     else if(x[2] < kCQ4Begin && x[2] > kCQ4End && rad2 < kCQ4SqRadius){
 174         b[0] = fQuadGradient*x[1];
 175         b[1] = fQuadGradient*x[0];
 176         b[2] = 0.;
 177     }
 178     else if(x[2] < kCD1Begin && x[2] > kCD1End && rad2 < kCD1SqRadius){
 179         b[1] = fDipoleField;
 180         b[2] = 0.;
 181         b[2] = 0.;
 182     }
 183     else if(x[2] < kCD2Begin && x[2] > kCD2End){
 184         if(((x[0]-kCD2XCentre1)*(x[0]-kCD2XCentre1)+(x[1]*x[1]))<kCD2SqRadius
 185            || ((x[0]-kCD2XCentre2)*(x[0]-kCD2XCentre2)+(x[1]*x[1]))<kCD2SqRadius){
 186           b[1] = -fDipoleField;
 187           b[2] = 0.;
 188           b[2] = 0.;
 189         }
 190     }
 191   }
 192
 193   // SIDE A **************************************************
 194   else{
 195     if(fCompensator && (x[2] > kACorr1Begin && x[2] < kACorr1End) && rad2 < kCCorr1SqRadius) {
 196       // Compensator magnet at z = 1075 m
 197         if (fFactor != 0.) {
 198             b[0] = fACorr1Field;
 199             b[1] = 0.;
 200             b[2] = 0.;
 201         }
 202         return;
 203     }
 204
 205     if(x[2] > kACorr2Begin && x[2] < kACorr2End && rad2 < kCCorr2SqRadius){
 206         if (fFactor != 0.) {
 207             b[0] = fACorr2Field;
 208             b[1] = 0.;
 209             b[2] = 0.;
 210         }
 211     }
 212     else if(x[2] > kAQ1Begin && x[2] < kAQ1End && rad2 < kAQ1SqRadius){
 213         // First quadrupole of inner triplet de-focussing in x-direction
 214         b[0] = -fQuadGradient*x[1];
 215         b[1] = -fQuadGradient*x[0];
 216         b[2] = 0.;
 217     }
 218     else if(x[2] > kAQ2Begin && x[2] < kAQ2End && rad2 < kAQ2SqRadius){
 219         b[0] = fQuadGradient*x[1];
 220         b[1] = fQuadGradient*x[0];
 221         b[2] = 0.;
 222     }
 223     else if(x[2] > kAQ3Begin && x[2] < kAQ3End && rad2 < kAQ3SqRadius){
 224         b[0] = fQuadGradient*x[1];
 225         b[1] = fQuadGradient*x[0];
 226         b[2] = 0.;
 227     }
 228     else if(x[2] > kAQ4Begin && x[2] < kAQ4End && rad2 < kAQ4SqRadius){
 229         b[0] = -fQuadGradient*x[1];
 230         b[1] = -fQuadGradient*x[0];
 231         b[2] = 0.;
 232     }
 233     else if(x[2] > kAD1Begin && x[2] < kAD1End && rad2 < kAD1SqRadius){
 234         b[0] = 0.;
 235         b[1] = -fDipoleField;
 236         b[2] = 0.;
 237     }
 238     else if(x[2] > kAD2Begin && x[2] < kAD2End){
 239         if(((x[0]-kAD2XCentre1)*(x[0]-kAD2XCentre1)+(x[1]*x[1])) < kAD2SqRadius
 240            || ((x[0]-kAD2XCentre2)*(x[0]-kAD2XCentre2)+(x[1]*x[1])) < kAD2SqRadius){
 241             b[1] = fDipoleField;
 242         }
 243     }
 244   }
 245 }
 246