HelicalHarmonicField Class Reference

Inheritance diagram for HelicalHarmonicField:

BLElementField

List of all members.


Detailed Description

HelicalHarmonicField represents one placement of a helicalharmonic magnet.

Public Member Functions

 HelicalHarmonicField (BLCoordinateTransform &_global2local, BLCMDhelicalharmonic *hh)
 constructor.
void addFieldValue (const G4double point[4], G4double field[6]) const
 addFieldValue() adds the field for this hh into field[]. point[] is in global coordinates.

Private Attributes

G4double radius
G4double halflength
int n
G4double b
G4double lambda
G4double phi0
BLCoordinateTransform global2local
G4RotationMatrix rotation

Constructor & Destructor Documentation

HelicalHarmonicField::HelicalHarmonicField ( BLCoordinateTransform _global2local,
BLCMDhelicalharmonic hh 
)

constructor.

References BLCMDhelicalharmonic::b, b, BLCoordinateTransform::getGlobal(), BLCoordinateTransform::getRotation(), global2local, halflength, BLCMDhelicalharmonic::lambda, lambda, BLCMDhelicalharmonic::length, BLCMDhelicalharmonic::n, n, BLCMDhelicalharmonic::phi0, phi0, BLCMDhelicalharmonic::radius, radius, rotation, and BLElementField::setGlobalPoint().

00307                                                                   :
00308                                         BLElementField(), rotation()
00309 {
00310         radius = hd->radius;
00311         halflength = hd->length/2.0;
00312         n = hd->n;
00313         b = hd->b;
00314         lambda = hd->lambda;
00315         phi0 = hd->phi0;
00316         global2local = _global2local;
00317         rotation = global2local.getRotation().inverse();
00318 
00319         // set global bounding box
00320         G4double local[4], global[4];
00321         local[3] = 0.0;
00322         for(int i=0; i<2; ++i) {
00323                 local[0] = (i==0 ? -1.0 : 1.0) * radius;
00324                 for(int j=0; j<2; ++j) {
00325                         local[1] = (j==0 ? -1.0 : 1.0) * radius;
00326                         for(int k=0; k<2; ++k) {
00327                                 local[2] = (k==0 ? -1.0 : 1.0) * halflength;
00328                                 global2local.getGlobal(local,global);
00329                                 setGlobalPoint(global);
00330                         }
00331                 }
00332         }
00333 }


Member Function Documentation

void HelicalHarmonicField::addFieldValue ( const G4double  point[4],
G4double  field[6] 
) const [virtual]

addFieldValue() adds the field for this hh into field[]. point[] is in global coordinates.

Implements BLElementField.

References b, B, BLCoordinateTransform::getLocal(), global2local, halflength, BLCoordinateTransform::isRotated(), lambda, n, phi0, radius, and rotation.

00336 {
00337   const char *CoDE="HelicalHarmonicField::addFieldValue";
00338   const int X=0,Y=1,Z=2;
00339   const int PHI=0,RHO=1;                      /* notation */
00340 
00341   G4ThreeVector Bcyl,B;
00342   G4ThreeVector global(point[X],point[Y],point[Z]);
00343   G4ThreeVector local;
00344 
00345   G4double phi,rho,cosphi,sinphi,k,kz,nkrho,psi,cospsi,sinpsi,Inm1,In,Inp1,b;
00346 
00347   global2local.getLocal(local,global);                  /*fetch local location in mm*/
00348   phi = atan2( local[Y], local[X]);
00349   rho = sqrt( local[X]*local[X] + local[Y]*local[Y] );  /*mm*/
00350   
00351   if( rho>radius || fabs(local[Z]) > halflength)  return;
00352 
00353         cosphi = cos(phi); 
00354         sinphi = sin(phi); 
00355         k = 2*pi/lambda; 
00356         kz = k*(local[Z]+halflength); 
00357         psi = n*(phi-kz+phi0); 
00358         cospsi = cos(psi); 
00359         sinpsi = sin(psi); 
00360         nkrho = n*k*rho;
00361         Inm1=gsl_sf_bessel_In(n-1,nkrho); 
00362         In=gsl_sf_bessel_In(n,nkrho); 
00363         Inp1=gsl_sf_bessel_In(n+1,nkrho); 
00364         b=HelicalHarmonicField::b*pow(2/(n*k),(n-1)); 
00365         Bcyl[PHI]=b*(Inm1-Inp1)*cospsi; 
00366         Bcyl[RHO]=b*(Inm1+Inp1)*sinpsi;
00367         B[Z]=-2*b*In*cospsi; 
00368         B[X]=Bcyl[RHO]*cosphi-Bcyl[PHI]*sinphi; 
00369         B[Y]=Bcyl[RHO]*sinphi+Bcyl[PHI]*cosphi;
00370 
00371   /* Rotation if applicable */
00372   
00373   if(global2local.isRotated())  B = rotation * B;
00374   
00375   field[0] += B[X];        /* update the field */
00376   field[1] += B[Y];
00377   field[2] += B[Z];
00378 }


Member Data Documentation

G4double HelicalHarmonicField::radius [private]

int HelicalHarmonicField::n [private]

G4double HelicalHarmonicField::b [private]

G4double HelicalHarmonicField::lambda [private]

G4double HelicalHarmonicField::phi0 [private]

G4RotationMatrix HelicalHarmonicField::rotation [private]


The documentation for this class was generated from the following file:
g4beamline