BLCMDpillbox Class Reference

Inheritance diagram for BLCMDpillbox:

BLElement BLCommand

List of all members.


Detailed Description

BLCMDpillbox implements a pillbox.

Each placement of a BLCMDpillbox creates a PillboxField that is linked into BLGlobalField. Each PillboxField sets its timeOffset via the tune particle (re-stepping through its TimingVol if necessary).

If fieldMapFile is non-null, it is read as a BLFieldMap and determines the peak field; both B and E are multiplied by maxGradient of the cavity.

Public Member Functions

 BLCMDpillbox ()
 Default constructor. Defines the command, etc.
virtual ~BLCMDpillbox ()
 Destructor.
 BLCMDpillbox (const BLCMDpillbox &r)
 Copy constructor.
BLElementclone ()
 clone()
G4String commandName ()
 commandName() returns "pillbox".
int command (BLArgumentVector &argv, BLArgumentMap &namedArgs)
 command() implements the pillbox command.
void defineNamedArgs ()
 defineNamedArgs() defines the named arguments of the command.
void argChanged ()
 argChanged() does internal computations after some arg changed
void construct (G4RotationMatrix *relativeRotation, G4ThreeVector relativePosition, G4LogicalVolume *parent, G4String parentName, G4RotationMatrix *parentRotation, G4ThreeVector parentPosition)
 construct() - construct the pillbox. Creates a new PillboxField and adds it to BLGlobalField.
G4double getLength ()
 getLength() returns the overallLength of the pillbox
G4double getWidth ()
 getWidth() returns the outer radius of the pillbox
G4double getHeight ()
 getHeight() returns the outer radius of the pillbox
G4bool isOK ()
 isOK() returns true only if all placed pillboxes of this type have had their timing set, and if the tuning (if used) converged.
void generatePoints (int npoints, std::vector< G4ThreeVector > &v)
 generatePoints() from BLElement
G4bool isOutside (G4ThreeVector &local, G4double tolerance)
 isOutside() from BLElement

Private Attributes

G4double maxGradient
G4String color
G4double frequency
G4double innerLength
G4double innerRadius
G4double pipeThick
G4double wallThick
G4double irisRadius
G4double collarRadialThick
G4double collarThick
G4double win1Thick
G4double win1OuterRadius
G4double win2Thick
G4String winMat
G4double phaseAcc
G4double skinDepth
G4double timingTolerance
G4double maxStep
G4String cavityMaterial
G4double timeOffset
G4double timeIncrement
G4String fieldMapFile
G4int kill
G4double overallLength
G4double outerRadius
G4double omega
G4double rFactor
G4double E2Bfactor
BLFieldMapfieldMap
std::vector< PillboxField * > pillboxField

Friends

class PillboxField


Constructor & Destructor Documentation

BLCMDpillbox::BLCMDpillbox (  ) 

Default constructor. Defines the command, etc.

References BLCMDTYPE_ELEMENT, cavityMaterial, collarRadialThick, collarThick, color, E2Bfactor, fieldMap, fieldMapFile, frequency, innerLength, innerRadius, irisRadius, kill, maxGradient, maxStep, omega, outerRadius, overallLength, phaseAcc, pipeThick, BLCommand::registerCommand(), rFactor, BLCommand::setDescription(), BLCommand::setSynopsis(), skinDepth, timeIncrement, timeOffset, timingTolerance, UNDETERMINED, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

Referenced by clone(), and command().

00189                            : BLElement(), pillboxField()
00190 {
00191         // register the commandName(), and its synopsis and description.
00192         registerCommand(BLCMDTYPE_ELEMENT);
00193         setSynopsis("Defines a pillbox RF cavity");
00194         setDescription("A Pillbox RF cavity is the basic RF element used to\n"
00195                 "construct a linac. The phaseAcc parameter sets the\n"
00196                 "phase of the tune particle at the center of the\n"
00197                 "cavity, and the timing offset of the cavity is determined\n"
00198                 "from that the first time that the Tune particle is tracked "
00199                 "through the cavity. Zero degrees is the rising zero-crossing "
00200                 "of the Ez field. If timeOffset is specified, it is used "
00201                 "rather than setting it from the Tune particle.\n\n"
00202                 "The Pipe, walls, and collars are always made of copper. "
00203                 "Pipe, wall, collar, win1, and win2 can be omitted by setting\n"
00204                 "their thickness to 0. Common usage is to set the collar "
00205                 "values so by placing multiple pillboxes sequentially the "
00206                 "collars form a beam pipe between them.\n\n"
00207                 "Note that section 4.4 of the User's "
00208                 "Guide has a dimensioned drawing of a pillbox.");
00209 
00210         // provide initial values for fields
00211         maxGradient = UNDETERMINED; // will be converted to MV/m
00212         color = "1.0,0.0,0.0";
00213         frequency = UNDETERMINED;
00214         innerLength = UNDETERMINED;
00215         innerRadius = UNDETERMINED;
00216         pipeThick = 3.0*mm;
00217         wallThick = 3.0*mm;
00218         irisRadius = 11.0*cm;
00219         collarRadialThick = 5.0*mm;
00220         collarThick = 2.5*mm;
00221         win1Thick = 0.200*mm;
00222         win1OuterRadius = 5.0*cm;
00223         win2Thick = 0.500*mm;
00224         winMat = "Be";
00225         phaseAcc = 40.0*deg;
00226         skinDepth = 0.002*mm;
00227         timingTolerance = 0.001*ns;
00228         maxStep = -1.0;
00229         cavityMaterial = "Vacuum";
00230         timeOffset = UNDETERMINED;
00231         timeIncrement = 0.0;
00232         fieldMapFile = "";
00233         kill = 0;
00234         // non-arguments:
00235         overallLength = UNDETERMINED;
00236         outerRadius = UNDETERMINED;
00237         omega = UNDETERMINED;
00238         rFactor = UNDETERMINED;
00239         E2Bfactor = UNDETERMINED;
00240         fieldMap = 0;
00241 }

virtual BLCMDpillbox::~BLCMDpillbox (  )  [inline, virtual]

Destructor.

00136 { }

BLCMDpillbox::BLCMDpillbox ( const BLCMDpillbox r  ) 

Copy constructor.

References cavityMaterial, collarRadialThick, collarThick, color, BLTune::copyTunableArg(), E2Bfactor, fieldMap, fieldMapFile, frequency, innerLength, innerRadius, irisRadius, kill, maxGradient, maxStep, omega, outerRadius, overallLength, phaseAcc, pipeThick, rFactor, skinDepth, timeIncrement, timeOffset, timingTolerance, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

00244                                                 : BLElement(r), pillboxField()
00245 {
00246         // copy fields one at a time (transfers default values from the
00247         // default object to this new object).
00248         BLTune::copyTunableArg(&maxGradient,&r.maxGradient);
00249         color = r.color;
00250         frequency = r.frequency;
00251         innerLength = r.innerLength;
00252         innerRadius = r.innerRadius;
00253         pipeThick = r.pipeThick;
00254         wallThick = r.wallThick;
00255         irisRadius = r.irisRadius;
00256         collarRadialThick = r.collarRadialThick;
00257         collarThick = r.collarThick;
00258         win1Thick = r.win1Thick;
00259         win1OuterRadius = r.win1OuterRadius;
00260         win2Thick = r.win2Thick;
00261         winMat = r.winMat;
00262         phaseAcc = r.phaseAcc;
00263         skinDepth = r.skinDepth;
00264         timingTolerance = r.timingTolerance;
00265         maxStep = r.maxStep;
00266         cavityMaterial = r.cavityMaterial;
00267         timeOffset = r.timeOffset;
00268         timeIncrement = r.timeIncrement;
00269         fieldMapFile = r.fieldMapFile;
00270         kill = r.kill;
00271 
00272         // non-arguments:
00273         overallLength = r.overallLength;
00274         outerRadius = r.outerRadius;
00275         omega = r.omega;
00276         rFactor = r.rFactor;
00277         E2Bfactor = r.E2Bfactor;
00278         fieldMap = r.fieldMap;
00279 }


Member Function Documentation

BLElement* BLCMDpillbox::clone (  )  [inline, virtual]

clone()

Implements BLElement.

References BLCMDpillbox().

00142 { return new BLCMDpillbox(*this); }

G4String BLCMDpillbox::commandName (  )  [inline, virtual]

commandName() returns "pillbox".

Implements BLCommand.

00145 { return "pillbox"; }

int BLCMDpillbox::command ( BLArgumentVector argv,
BLArgumentMap namedArgs 
) [virtual]

command() implements the pillbox command.

Implements BLCommand.

References argChanged(), BLCMDpillbox(), cavityMaterial, E, E2Bfactor, fieldMap, fieldMapFile, frequency, g4bl_exit(), BLCommand::getMaterial(), BLCommand::handleNamedArgs(), innerRadius, maxGradient, BLCommand::print(), BLCommand::printError(), BLFieldMap::readFile(), and BLElement::setName().

00282 {
00283         // If the name of the cavity is "testUnits", perform the test and exit.
00284         if(argv[0] == "testUnits") {
00285                 // set the test case values, from Al Moretti's SuperFish example
00286                 maxGradient = 1.0 * megavolt/meter;
00287                 frequency = 0.956188/ns;
00288                 double Hcorrect = 1545.0; // A/m
00289                 double BcorrectTesla = 4.0*pi*1E-7 * Hcorrect; // MKS
00290                 argChanged();
00291                 // the first maximum of J1() has a value of 0.5819
00292                 double Bmax = maxGradient*E2Bfactor*0.5819;
00293                 printf("pillobx testUnits:\n           "         
00294                         "Freq(GHz)  Emax(MV/m)  Bmax(Tesla)  Radius(mm)\n");
00295                 printf("Values:  %10.4f %10.4f %12.4e %11.4f\n",
00296                         frequency/(1e9*hertz), maxGradient/(megavolt/meter),
00297                         Bmax/tesla, innerRadius/mm);
00298                 // these values already are in the correct units:
00299                 printf("Correct: %10.4f %10.4f %12.4e %8.1f\n",
00300                         0.956188,1.0,BcorrectTesla, 120.0);
00301                 printf("pillbox testUnits: exiting\n");
00302                 extern void g4bl_exit(int);
00303                 g4bl_exit(0);
00304         }
00305 
00306         if(argv.size() != 1) {
00307                 printError("pillbox: Invalid command, must have name");
00308                 return -1;
00309         }
00310 
00311         if(argv[0] == "default") {
00312                 return defaultPillbox.handleNamedArgs(namedArgs);
00313         }
00314 
00315         BLCMDpillbox *p = new BLCMDpillbox(defaultPillbox);
00316         p->setName(argv[0]);
00317         int retval = p->handleNamedArgs(namedArgs);
00318 
00319         if(p->fieldMapFile != "") {
00320                 p->fieldMap = new BLFieldMap();
00321                 p->fieldMap->readFile(p->fieldMapFile);
00322         }
00323 
00324         // check material exists
00325         if(p->cavityMaterial.size() > 0) getMaterial(p->cavityMaterial,false);
00326 
00327         p->print(argv[0]);
00328 
00329         return retval;
00330 }

void BLCMDpillbox::defineNamedArgs (  )  [virtual]

defineNamedArgs() defines the named arguments of the command.

Reimplemented from BLCommand.

References BLCommand::argDouble(), BLCommand::argInt(), BLCommand::argString(), BLCommand::argTunable(), cavityMaterial, collarRadialThick, collarThick, color, fieldMapFile, frequency, innerLength, innerRadius, irisRadius, kill, maxGradient, maxStep, phaseAcc, pipeThick, skinDepth, timeIncrement, timeOffset, timingTolerance, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

00333 {
00334         argTunable(maxGradient,"maxGradient","The peak gradient of the cavity (MV/m)",megavolt/meter);
00335         argString(color,"color","The color of the cavity");
00336         argDouble(frequency,"frequency","The frequency of the cavity (GHz)",1/ns,"",false);
00337         argDouble(innerLength,"innerLength","The inside length of the cavity (mm)",mm,"",false);
00338         argDouble(innerRadius,"innerRadius","The inside radius of the cavity (mm)",mm,"",false);
00339         argDouble(pipeThick,"pipeThick","The thickness of the pipe wall (mm)",mm,"",false);
00340         argDouble(wallThick,"wallThick","The thickness of the cavity walls (mm)",mm,"",false);
00341         argDouble(irisRadius,"irisRadius","The radius of the iris (mm)",mm,"",false);
00342         argDouble(collarRadialThick,"collarRadialThick","The radial thickness of the collar (mm)",mm,"",false);
00343         argDouble(collarThick,"collarThick","The thickness of the collar along z(mm)",mm,"",false);
00344         argDouble(win1Thick,"win1Thick","The thickness of the central portion of the\n"
00345                 "windows; zero for no window (mm)",mm,"",false);
00346         argDouble(win1OuterRadius,"win1OuterRadius","The radius of the central portion of\n"
00347                 "the windows (mm)",mm,"",false);
00348         argDouble(win2Thick,"win2Thick","The thickness of the outer portion of the\n"
00349                 "windows; zero for no window (mm)",mm,"",false);
00350         argString(winMat,"winMat","The material of the windows");
00351         argDouble(phaseAcc,"phaseAcc","The reference phase of the cavity (degrees)",deg);
00352         argDouble(skinDepth,"skinDepth","The skin depth (mm)",mm,"",false);
00353         argDouble(timingTolerance,"timingTolerance","Tolerance for timing tuning (ns)",ns);
00354         argDouble(maxStep,"maxStep","The maximum stepsize in the element (mm).",mm);
00355         argString(cavityMaterial,"cavityMaterial","Material of cavity volume (Vacuum).");
00356         argDouble(timeOffset,"timeOffset","Time offset for cavity (default: tuned by tune particle) (ns).",ns);
00357         argDouble(timeIncrement,"timeIncrement","Increment to timeOffset, applied AFTER tuning. (ns).",ns);
00358         argString(fieldMapFile,"fieldMapFile","Filename for BLFieldMap (pillbox if null).",false);
00359         argInt(kill,"kill","Set nonzero to kill tracks that hit the pipe, walls, or collars (0).");
00360 }

void BLCMDpillbox::argChanged (  )  [virtual]

argChanged() does internal computations after some arg changed

Reimplemented from BLCommand.

References collarThick, E2Bfactor, frequency, BLParam::getDouble(), innerLength, innerRadius, maxStep, omega, outerRadius, overallLength, Param, pipeThick, BLCommand::printError(), rFactor, UNDETERMINED, wallThick, win1Thick, and win2Thick.

Referenced by command().

00363 {
00364         // compute non-argument variables
00365         omega = 2.0 * pi * frequency;
00366         if(innerRadius == UNDETERMINED) {
00367                 if(frequency > 0.0) {
00368                         // the first zero of J0(r) is at r=2.405
00369                         innerRadius =  2.405 / (omega/c_light);
00370                 } else {
00371                         printError("pillbox: innerRadius is undefined");
00372                 }
00373         }
00374         double wall=collarThick;
00375         if(wallThick > wall) wall = wallThick;
00376         if(win1Thick > wall) wall = win1Thick;
00377         if(win2Thick > wall) wall = win2Thick;
00378         overallLength = innerLength + 2.0 * wall;
00379         outerRadius = innerRadius + pipeThick;
00380         rFactor = omega/c_light;
00381         // Previous versions had the incorrect value from the BeamTools:
00382         //     E2Bfactor = 1.0e-9 / c_light;   // MeV/m, MKS units
00383         // C.f. Eq. 8.77 in J.D.Jackson, _Classical_Electrodynamics_,
00384         //       first edition (1962). His value is 1/c in vacuum.
00385         E2Bfactor = 1.0/c_light;   // Geant4 units
00386         if(maxStep < 0.0) maxStep = Param.getDouble("maxStep");
00387 }

void BLCMDpillbox::construct ( G4RotationMatrix *  relativeRotation,
G4ThreeVector  relativePosition,
G4LogicalVolume *  parent,
G4String  parentName,
G4RotationMatrix *  parentRotation,
G4ThreeVector  parentPosition 
) [virtual]

construct() - construct the pillbox. Creates a new PillboxField and adds it to BLGlobalField.

Implements BLElement.

References BLGlobalField::addElementField(), cavityMaterial, collarRadialThick, collarThick, BLCommand::getMaterial(), BLElement::getName(), BLGlobalField::getObject(), BLManager::getObject(), BLCoordinateTransform::getPosition(), PillboxField::global2local, innerLength, innerRadius, irisRadius, BLCoordinateTransform::isRotated(), kill, maxStep, PillboxField::pillbox, pillboxField, PillboxField, pipeThick, BLManager::registerTuneParticleStep(), PillboxField::timingPV, wallThick, win1OuterRadius, win1Thick, win2Thick, and winMat.

00395 {
00396         G4String thisname = parentName+getName();
00397 
00398         G4Material *cu = getMaterial("Cu");
00399         G4Material *wm = getMaterial(winMat);
00400         G4Material *cavityMat = getMaterial(cavityMaterial);
00401         G4UserLimits *userLimits = new G4UserLimits(maxStep);
00402 
00403         // find enclosing Tubs
00404         G4double halflen = innerLength/2+wallThick;
00405         halflen = std::max(halflen,innerLength/2+collarThick);
00406         halflen = std::max(halflen,innerLength/2+win1Thick);
00407         halflen = std::max(halflen,innerLength/2+win2Thick);
00408         G4Tubs *tubs = new G4Tubs(getName(),0.0,innerRadius+pipeThick,
00409                                         halflen,0.0,2.0*pi);
00410         G4LogicalVolume *pillbox = new G4LogicalVolume(tubs,cavityMat,
00411                                         parentName+getName(),0,0,userLimits);
00412 
00413         G4ThreeVector loc;
00414         G4VPhysicalVolume *pv;
00415 
00416         // The outer pipe
00417         G4LogicalVolume *pipe=0;
00418         if(pipeThick > 0.001*mm) {
00419                 tubs = new G4Tubs(getName()+"pipe",innerRadius,
00420                                         innerRadius+pipeThick,
00421                                         innerLength/2+wallThick,0.0,2.0*pi);
00422                 pipe = new G4LogicalVolume(tubs,cu,
00423                                         getName()+"pipeLog",0,0,userLimits);
00424                 pv = new G4PVPlacement(0,loc,pipe,parentName+getName()+"Pipe",
00425                                         pillbox,false,0);
00426                 if(kill) BLManager::getObject()->
00427                         registerSteppingAction(pv,new BLKillTrack(thisname));
00428         }
00429 
00430         // the wall (us and ds) - radially: collarOuterRadius to innerRadius
00431         G4LogicalVolume *wall=0;
00432         if(wallThick > 0.001*mm) {
00433                 tubs = new G4Tubs(getName()+"wall",irisRadius+collarRadialThick,
00434                                         innerRadius,wallThick/2,0.0,2.0*pi);
00435                 wall = new G4LogicalVolume(tubs,cu,
00436                                         getName()+"wallLog",0,0,userLimits);
00437                 G4double wallCenterZ = innerLength/2.0 + wallThick/2.0;
00438                 loc.setZ(-wallCenterZ);
00439                 pv = new G4PVPlacement(0,loc,wall,parentName+getName()+"UsWall",
00440                                         pillbox,false,0);
00441                 if(kill) BLManager::getObject()->
00442                         registerSteppingAction(pv,new BLKillTrack(thisname));
00443                 loc.setZ(+wallCenterZ);
00444                 pv = new G4PVPlacement(0,loc,wall,parentName+getName()+"DsWall",
00445                                         pillbox,false,0);
00446                 if(kill) BLManager::getObject()->
00447                         registerSteppingAction(pv,new BLKillTrack(thisname));
00448         }
00449         
00450         // the collar (us and ds) - includes wall thickness (i.e. the inner
00451         // edge of the collar is the same as the inner edge of the wall)
00452         G4LogicalVolume *collar=0;
00453         if(collarThick > 0.001*mm) {
00454                 tubs = new G4Tubs(getName()+"collar",irisRadius,
00455                                         irisRadius+collarRadialThick,
00456                                         collarThick/2,0.0,2.0*pi);
00457                 collar = new G4LogicalVolume(tubs,cu,
00458                                         getName()+"collarLog",0,0,userLimits);
00459                 G4double collarCenterZ = innerLength/2.0 + collarThick/2.0;
00460                 loc.setZ(-collarCenterZ);
00461                 pv = new G4PVPlacement(0,loc,collar,parentName+getName()+"UsCollar",
00462                                         pillbox,false,0);
00463                 if(kill) BLManager::getObject()->
00464                         registerSteppingAction(pv,new BLKillTrack(thisname));
00465                 loc.setZ(+collarCenterZ);
00466                 pv = new G4PVPlacement(0,loc,collar,parentName+getName()+"DsCollar",
00467                                         pillbox,false,0);
00468                 if(kill) BLManager::getObject()->
00469                         registerSteppingAction(pv,new BLKillTrack(thisname));
00470         }
00471         
00472         // win1 (the inner window)
00473         G4LogicalVolume *win1=0;
00474         if(win1Thick > 0.001*mm) {
00475                 tubs = new G4Tubs(getName()+"win1",0.0,win1OuterRadius,
00476                                         win1Thick/2,0.0,2.0*pi);
00477                 win1 = new G4LogicalVolume(tubs,wm,
00478                                         getName()+"win1Log",0,0,userLimits);
00479                 G4double win1CenterZ = innerLength/2.0 + win1Thick/2.0;
00480                 loc.setZ(-win1CenterZ);
00481                 new G4PVPlacement(0,loc,win1,parentName+getName()+"UsWin1",
00482                                         pillbox,false,0);
00483                 loc.setZ(+win1CenterZ);
00484                 new G4PVPlacement(0,loc,win1,parentName+getName()+"DsWin1",
00485                                         pillbox,false,0);
00486         }
00487 
00488         // win2 (the outer window)
00489         G4LogicalVolume *win2=0;
00490         if(win2Thick > 0.001*mm) {
00491                 tubs = new G4Tubs(getName()+"win2",win1OuterRadius,irisRadius,
00492                                         win2Thick/2,0.0,2.0*pi);
00493                 win2 = new G4LogicalVolume(tubs,wm,
00494                                         getName()+"win2Log",0,0,userLimits);
00495                 G4double win2CenterZ = innerLength/2.0 + win2Thick/2.0;
00496                 loc.setZ(-win2CenterZ);
00497                 new G4PVPlacement(0,loc,win2,parentName+getName()+"UsWin2",
00498                                         pillbox,false,0);
00499                 loc.setZ(+win2CenterZ);
00500                 new G4PVPlacement(0,loc,win2,parentName+getName()+"DsWin2",
00501                                         pillbox,false,0);
00502         }
00503 
00504         // half of the interior
00505         tubs = new G4Tubs(getName()+"half",0.0,innerRadius,
00506                                         innerLength/4,0.0,2.0*pi);
00507         G4LogicalVolume *half = new G4LogicalVolume(tubs,cavityMat,getName()+"halfLog",
00508                                         0,0,userLimits);
00509         G4double halfCenterZ = innerLength/4.0;
00510         loc.setZ(+halfCenterZ);
00511         new G4PVPlacement(0,loc,half,parentName+getName()+"DsHalf",
00512                                                         pillbox,false,0);
00513         loc.setZ(-halfCenterZ);
00514         G4VPhysicalVolume *timingPV = new G4PVPlacement(0,loc,half,
00515                                 parentName+getName()+"TimingVol",pillbox,false,0);
00516 
00517 #ifdef G4VIS_USE
00518         G4VisAttributes *c1 = new G4VisAttributes(true,G4Color(0.5,0.0,0.0));
00519         G4VisAttributes *c2 = new G4VisAttributes(true,G4Color(0.7,0.0,0.0));
00520         G4VisAttributes *c3 = new G4VisAttributes(true,G4Color(1.0,0.0,0.0));
00521         const G4VisAttributes *invisible = &G4VisAttributes::GetInvisible();
00522         if(pipe) pipe->SetVisAttributes(c1);
00523         if(wall) wall->SetVisAttributes(c1);
00524         if(collar) collar->SetVisAttributes(c1);
00525         if(win1) win1->SetVisAttributes(c3);
00526         if(win2) win2->SetVisAttributes(c2);
00527         half->SetVisAttributes(invisible);
00528         pillbox->SetVisAttributes(invisible);
00529 #endif
00530 
00531         // geant4 rotation convention is backwards from g4beamline
00532         G4RotationMatrix *g4rot = 0;
00533         if(relativeRotation)
00534                 g4rot = new G4RotationMatrix(relativeRotation->inverse());
00535 
00536         // place the pillbox into the parent
00537         new G4PVPlacement(g4rot,relativePosition,pillbox,
00538                                 parentName+getName(),parent,false,0);
00539 
00540         // get globalRotation and globalPosition
00541         G4RotationMatrix *globalRotation = 0;
00542         if(relativeRotation && parentRotation) {
00543                 globalRotation = 
00544                     new G4RotationMatrix(*parentRotation * *relativeRotation);
00545         } else if(relativeRotation) {
00546                 globalRotation = relativeRotation;
00547         } else if(parentRotation) {
00548                 globalRotation = parentRotation;
00549         }
00550         G4ThreeVector globalPosition(relativePosition + parentPosition);
00551         if(parentRotation)
00552                 globalPosition = *parentRotation * relativePosition +
00553                                 parentPosition;
00554 
00555         BLCoordinateTransform global2local(globalRotation,globalPosition);
00556 
00557         // add this pillbox field to the GlobalField, and to BLManager
00558         G4double dz = innerLength/2.0 + wallThick;
00559         // if rotated, make an overestimate of the field occupancy along z
00560         if(global2local.isRotated())
00561                 dz += innerRadius;
00562         PillboxField *pf = 
00563                 new PillboxField(thisname,global2local,this,timingPV);
00564         BLGlobalField::getObject()->addElementField(pf);
00565         pillboxField.push_back(pf);
00566         BLManager::getObject()->registerTuneParticleStep(timingPV,pf);
00567 
00568         printf("BLCMDpillbox::construct %s parent=%s relZ=%.1f globZ=%.1f\n"
00569                         "\tzmin=%.1f zmax=%.1f\n",
00570                 thisname.c_str(),parentName.c_str(),relativePosition[2],
00571                 global2local.getPosition()[2], global2local.getPosition()[2]-dz,
00572                 global2local.getPosition()[2]+dz);
00573 }

G4double BLCMDpillbox::getLength (  )  [inline, virtual]

getLength() returns the overallLength of the pillbox

Implements BLElement.

References overallLength.

00166 { return overallLength; }

G4double BLCMDpillbox::getWidth (  )  [inline, virtual]

getWidth() returns the outer radius of the pillbox

Implements BLElement.

References outerRadius.

00169 { return outerRadius*2.0; }

G4double BLCMDpillbox::getHeight (  )  [inline, virtual]

getHeight() returns the outer radius of the pillbox

Implements BLElement.

References outerRadius.

00172 { return outerRadius*2.0; }

G4bool BLCMDpillbox::isOK (  )  [virtual]

isOK() returns true only if all placed pillboxes of this type have had their timing set, and if the tuning (if used) converged.

Implements BLElement.

References pillboxField, and PillboxField::TIMING_COMPLETE.

00576 {
00577         G4bool retval = true;
00578 
00579         // verify all placed pillboxField-s have set their timeOffset 
00580         std::vector<PillboxField *>::iterator i;
00581         for(i=pillboxField.begin(); i!=pillboxField.end(); ++i) {
00582                 if((*i)->state != PillboxField::TIMING_COMPLETE) {
00583                         printf("*** BLCMDpillbox %s has not determined timing\n",
00584                                 (*i)->getName().c_str());
00585                         retval = false;
00586                 }
00587         }
00588 
00589         return retval;
00590 }

void BLCMDpillbox::generatePoints ( int  npoints,
std::vector< G4ThreeVector > &  v 
) [virtual]

generatePoints() from BLElement

Implements BLElement.

References BLElement::generateTubs(), outerRadius, and overallLength.

00784 {
00785         generateTubs(npoints, 0.0, outerRadius, 0.0, 360.0*deg,
00786                         overallLength, v);
00787 }

G4bool BLCMDpillbox::isOutside ( G4ThreeVector &  local,
G4double  tolerance 
) [virtual]

isOutside() from BLElement

Implements BLElement.

References outerRadius, and overallLength.

00790 {
00791         G4double r = sqrt(local[0]*local[0]+local[1]*local[1]);
00792         return r > outerRadius-tolerance ||
00793                 fabs(local[2]) > overallLength/2.0-tolerance;
00794 }


Friends And Related Function Documentation

friend class PillboxField [friend]

Referenced by construct().


Member Data Documentation

G4double BLCMDpillbox::maxGradient [private]

G4String BLCMDpillbox::color [private]

Referenced by BLCMDpillbox(), and defineNamedArgs().

G4double BLCMDpillbox::frequency [private]

G4double BLCMDpillbox::innerLength [private]

G4double BLCMDpillbox::innerRadius [private]

G4double BLCMDpillbox::pipeThick [private]

G4double BLCMDpillbox::wallThick [private]

G4double BLCMDpillbox::irisRadius [private]

G4double BLCMDpillbox::collarRadialThick [private]

G4double BLCMDpillbox::collarThick [private]

G4double BLCMDpillbox::win1Thick [private]

G4double BLCMDpillbox::win1OuterRadius [private]

G4double BLCMDpillbox::win2Thick [private]

G4String BLCMDpillbox::winMat [private]

G4double BLCMDpillbox::phaseAcc [private]

G4double BLCMDpillbox::skinDepth [private]

G4double BLCMDpillbox::timingTolerance [private]

G4double BLCMDpillbox::maxStep [private]

G4String BLCMDpillbox::cavityMaterial [private]

G4double BLCMDpillbox::timeOffset [private]

G4double BLCMDpillbox::timeIncrement [private]

G4String BLCMDpillbox::fieldMapFile [private]

G4int BLCMDpillbox::kill [private]

G4double BLCMDpillbox::overallLength [private]

G4double BLCMDpillbox::outerRadius [private]

G4double BLCMDpillbox::omega [private]

G4double BLCMDpillbox::rFactor [private]

G4double BLCMDpillbox::E2Bfactor [private]

std::vector<PillboxField *> BLCMDpillbox::pillboxField [private]

Referenced by construct(), and isOK().


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