frequency=0.0 is accepted and yields a pillbox with a constant Ez = maxGradient (useful to verify units are correct).
Public Member Functions | |
PillboxField (G4String &_name, BLCoordinateTransform &_global2local, BLCMDpillbox *_pillbox, G4VPhysicalVolume *_timingPV) | |
constructor. _zmin/max are global coordinates. | |
G4String | getName () |
getName() returns the name of this placement. | |
void | addFieldValue (const G4double point[4], G4double field[6]) const |
addFieldValue() adds the field for this pillbox into field[]. point[] is in global coordinates. | |
void | UserSteppingAction (const G4Step *step) |
doStep() handles a step of the tune particle. | |
Private Types | |
enum | State { TIMING_UNKNOWN, SETTING_TIMING, TIMING_COMPLETE } |
Private Attributes | |
G4String | name |
BLCMDpillbox * | pillbox |
G4VPhysicalVolume * | timingPV |
BLCoordinateTransform | global2local |
G4RotationMatrix * | rotation |
G4double | timeOffset |
G4Track | saveTrack |
G4int | timeCount |
G4bool | validSaveTrack |
State | state |
BLFieldMap * | fieldMap |
Friends | |
class | BLCMDpillbox |
enum PillboxField::State [private] |
PillboxField::PillboxField | ( | G4String & | _name, | |
BLCoordinateTransform & | _global2local, | |||
BLCMDpillbox * | _pillbox, | |||
G4VPhysicalVolume * | _timingPV | |||
) |
constructor. _zmin/max are global coordinates.
References BLCMDpillbox::fieldMap, fieldMap, BLCMDpillbox::frequency, BLCoordinateTransform::getGlobal(), BLCoordinateTransform::getRotation(), global2local, BLCMDpillbox::innerLength, BLCMDpillbox::innerRadius, BLCoordinateTransform::isRotated(), isUndetermined(), BLCMDpillbox::maxGradient, name, pillbox, rotation, BLElementField::setGlobalPoint(), state, timeCount, BLCMDpillbox::timeOffset, timeOffset, TIMING_COMPLETE, TIMING_UNKNOWN, timingPV, validSaveTrack, and BLCMDpillbox::wallThick.
00595 : 00596 BLElementField(), BLManager::SteppingAction() 00597 { 00598 name = _name; 00599 pillbox = _pillbox; 00600 timingPV = _timingPV; 00601 global2local = _global2local; 00602 rotation = 0; 00603 timeOffset = pillbox->timeOffset; 00604 timeCount = 0; 00605 validSaveTrack = false; 00606 state = TIMING_UNKNOWN; 00607 fieldMap = pillbox->fieldMap; 00608 00609 if(global2local.isRotated()) { 00610 rotation = new G4RotationMatrix(global2local.getRotation()); 00611 rotation->invert(); 00612 } 00613 00614 if(pillbox->frequency == 0.0 || pillbox->maxGradient == 0.0) 00615 timeOffset = 0.0; 00616 00617 // don't tune timeOffset if it is known 00618 if(!isUndetermined(timeOffset)) 00619 state = TIMING_COMPLETE; 00620 00621 // set global bounding box 00622 G4double local[4], global[4]; 00623 local[3] = 0.0; 00624 G4double dz = pillbox->innerLength/2.0 + pillbox->wallThick; 00625 for(int i=0; i<2; ++i) { 00626 local[0] = (i==0 ? -1.0 : 1.0) * pillbox->innerRadius; 00627 for(int j=0; j<2; ++j) { 00628 local[1] = (j==0 ? -1.0 : 1.0) * pillbox->innerRadius; 00629 for(int k=0; k<2; ++k) { 00630 local[2] = (k==0 ? -1.0 : 1.0) * dz; 00631 global2local.getGlobal(local,global); 00632 setGlobalPoint(global); 00633 } 00634 } 00635 } 00636 }
G4String PillboxField::getName | ( | ) | [inline] |
getName() returns the name of this placement.
References name.
Referenced by UserSteppingAction().
00080 { return name; }
void PillboxField::addFieldValue | ( | const G4double | point[4], | |
G4double | field[6] | |||
) | const [virtual] |
addFieldValue() adds the field for this pillbox into field[]. point[] is in global coordinates.
Implements BLElementField.
References B, E, BLCMDpillbox::E2Bfactor, fieldMap, BLCMDpillbox::frequency, BLFieldMap::getFieldValue(), BLCoordinateTransform::getLocal(), global2local, BLCMDpillbox::innerLength, BLCMDpillbox::innerRadius, BLCMDpillbox::maxGradient, BLCMDpillbox::omega, pillbox, BLCMDpillbox::rFactor, rotation, BLCMDpillbox::skinDepth, timeOffset, and BLCMDpillbox::wallThick.
00639 { 00640 G4double local[4]; 00641 global2local.getLocal(local,point); 00642 00643 if(fieldMap) { 00644 G4double myField[6]; 00645 fieldMap->getFieldValue(local,myField,pillbox->maxGradient, 00646 pillbox->maxGradient); 00647 // frequency == 0.0 means a d.c. electric field 00648 G4double timeE=1.0, timeB=0.0; 00649 if(pillbox->frequency != 0) { 00650 double tArg = pillbox->omega * (point[3] - timeOffset); 00651 timeE = cos(tArg); 00652 timeB = sin(tArg); 00653 } 00654 if(rotation) { 00655 G4ThreeVector B(myField[0],myField[1],myField[2]); 00656 G4ThreeVector E(myField[3],myField[4],myField[5]); 00657 B = *rotation * B; 00658 E = *rotation * E; 00659 field[0] += B[0] * timeB; 00660 field[1] += B[1] * timeB; 00661 field[2] += B[2] * timeB; 00662 field[3] += E[0] * timeE; 00663 field[4] += E[1] * timeE; 00664 field[5] += E[2] * timeE; 00665 } else { 00666 field[0] += myField[0] * timeB; 00667 field[1] += myField[1] * timeB; 00668 field[2] += myField[2] * timeB; 00669 field[3] += myField[3] * timeE; 00670 field[4] += myField[4] * timeE; 00671 field[5] += myField[5] * timeE; 00672 } 00673 return; 00674 } 00675 00676 G4double x = local[0]; 00677 G4double y = local[1]; 00678 G4double z = local[2]; 00679 G4double r = sqrt(x*x + y*y); 00680 G4double dz = pillbox->innerLength/2.0 + pillbox->wallThick; 00681 if(z < -dz || z > dz || r > pillbox->innerRadius) return; 00682 00683 // frequency == 0.0 means a d.c. electric field 00684 if(pillbox->frequency == 0.0) { 00685 // zmin and zmax include the walls, so re-test z position 00686 if(z >= -pillbox->innerLength/2.0 && 00687 z <= pillbox->innerLength/2.0) 00688 field[5] += pillbox->maxGradient; 00689 return; 00690 } 00691 00692 // compute normal RF pillbox field 00693 double arg = pillbox->rFactor*r; 00694 if(fabs(arg) < 1.0e-20) arg = 0.0; 00695 double j0 = gsl_sf_bessel_J0(arg); 00696 double j1 = gsl_sf_bessel_J1(arg); 00697 double tArg = pillbox->omega * (point[3] - timeOffset); 00698 double ez = pillbox->maxGradient * j0 * cos(tArg); 00699 double bphi = -pillbox->maxGradient*pillbox->E2Bfactor*j1*sin(tArg); 00700 // handle skin depth if in wall or window 00701 double f = 1.0; 00702 if(z > pillbox->innerLength/2.0) 00703 f=exp(-(z-pillbox->innerLength/2.0)/pillbox->skinDepth); 00704 else if(z < -pillbox->innerLength/2.0) 00705 f=exp(-(-z-pillbox->innerLength/2.0)/pillbox->skinDepth); 00706 double phi = atan2(point[1],point[0]); 00707 if(rotation) { 00708 G4ThreeVector B(-bphi*sin(phi)*f, bphi*cos(phi)*f, 0.0); 00709 G4ThreeVector E(0.0, 0.0, ez*f); 00710 B = *rotation * B; 00711 E = *rotation * E; 00712 field[0] += B[0]; 00713 field[1] += B[1]; 00714 field[2] += B[2]; 00715 field[3] += E[0]; 00716 field[4] += E[1]; 00717 field[5] += E[2]; 00718 } else { 00719 field[0] += -bphi * sin(phi) * f; 00720 field[1] += bphi * cos(phi) * f; 00721 field[5] += ez * f; 00722 } 00723 }
void PillboxField::UserSteppingAction | ( | const G4Step * | step | ) | [virtual] |
doStep() handles a step of the tune particle.
Implements BLManager::SteppingAction.
References BLCMDpillbox::frequency, getName(), BLManager::getObject(), BLManager::getSteppingManager(), BLCMDpillbox::innerLength, ITERATION_LIMIT, BLCMDpillbox::maxGradient, BLCMDpillbox::omega, BLCMDpillbox::phaseAcc, pillbox, saveTrack, SETTING_TIMING, state, timeCount, BLCMDpillbox::timeIncrement, timeOffset, TIMING_COMPLETE, timingPV, BLCMDpillbox::timingTolerance, and validSaveTrack.
00727 { 00728 // do nothing if d.c. or no field or timing is complete 00729 if(pillbox->frequency <= 0.0 || pillbox->maxGradient == 0.0) 00730 state = TIMING_COMPLETE; 00731 if(state == TIMING_COMPLETE) return; 00732 00733 G4Track *track = step->GetTrack(); 00734 G4StepPoint *prePoint = step->GetPreStepPoint(); 00735 if(!prePoint) return; 00736 G4VPhysicalVolume *prePV = prePoint->GetPhysicalVolume(); 00737 G4StepPoint *postPoint = step->GetPostStepPoint(); 00738 if(!postPoint) return; 00739 G4VPhysicalVolume *postPV = postPoint->GetPhysicalVolume(); 00740 G4SteppingManager *steppingMgr = 00741 BLManager::getObject()->getSteppingManager(); 00742 00743 if(prePV == postPV) return; // neither entering nor leaving 00744 00745 if(postPV == timingPV) { // entering timingPV 00746 // save track for setting timeOffset 00747 saveTrack.CopyTrackInfo(*step->GetTrack()); 00748 saveTrack.SetUserInformation(0); 00749 validSaveTrack = true; 00750 state = SETTING_TIMING; 00751 // estimate timingOffset by neglecting acceleration 00752 double v = saveTrack.GetVelocity(); 00753 double dist = pillbox->innerLength/2.0; 00754 double arrival = saveTrack.GetGlobalTime() + dist/v; 00755 // pi/2 because E is a cos() and we want phaseAcc=0 00756 // to be the rising zero-crossing 00757 timeOffset = arrival - (pillbox->phaseAcc-pi/2.0)/pillbox->omega; 00758 } else if(prePV == timingPV) { // leaving timingPV 00759 if(++timeCount > ITERATION_LIMIT) 00760 G4Exception("pillbox","Iteration Limit",FatalException,""); 00761 // check timing offset 00762 double t = track->GetGlobalTime(); 00763 double phase = pillbox->omega * (t-timeOffset) + pi/2.0; 00764 double dt = (phase-pillbox->phaseAcc)/pillbox->omega; 00765 timeOffset = timeOffset + dt; 00766 if(fabs(dt) <= pillbox->timingTolerance) { 00767 // Success! 00768 printf("pillbox %s: Time OK timeOffset=%.4f ns, incremented to " 00769 "%.4f ns\n",getName().c_str(),timeOffset, 00770 timeOffset+pillbox->timeIncrement); 00771 timeOffset += pillbox->timeIncrement; 00772 state = TIMING_COMPLETE; 00773 return; 00774 } 00775 // restore saved data (i.e. jump back to when the track 00776 // enteered TimingVol) 00777 if(!validSaveTrack) 00778 G4Exception("pillbox","Invalid Step",FatalException,""); 00779 steppingMgr->GetfSecondary()->push_back(new G4Track(saveTrack)); 00780 track->SetTrackStatus(fStopAndKill); 00781 } 00782 }
friend class BLCMDpillbox [friend] |
G4String PillboxField::name [private] |
Referenced by getName(), and PillboxField().
BLCMDpillbox* PillboxField::pillbox [private] |
Referenced by addFieldValue(), BLCMDpillbox::construct(), PillboxField(), and UserSteppingAction().
G4VPhysicalVolume* PillboxField::timingPV [private] |
Referenced by BLCMDpillbox::construct(), PillboxField(), and UserSteppingAction().
Referenced by addFieldValue(), BLCMDpillbox::construct(), and PillboxField().
G4RotationMatrix* PillboxField::rotation [private] |
Referenced by addFieldValue(), and PillboxField().
G4double PillboxField::timeOffset [private] |
Referenced by addFieldValue(), PillboxField(), and UserSteppingAction().
G4Track PillboxField::saveTrack [private] |
Referenced by UserSteppingAction().
G4int PillboxField::timeCount [private] |
Referenced by PillboxField(), and UserSteppingAction().
G4bool PillboxField::validSaveTrack [private] |
Referenced by PillboxField(), and UserSteppingAction().
State PillboxField::state [private] |
Referenced by PillboxField(), and UserSteppingAction().
BLFieldMap* PillboxField::fieldMap [private] |
Referenced by addFieldValue(), and PillboxField().