#include <G4ParameterisationPolyhedra.hh>
Inheritance diagram for G4ParameterisationPolyhedraZ:
Public Member Functions | |
G4ParameterisationPolyhedraZ (EAxis axis, G4int nCopies, G4double offset, G4double step, G4VSolid *motherSolid, DivisionType divType) | |
~G4ParameterisationPolyhedraZ () | |
void | CheckParametersValidity () |
G4double | GetMaxParameter () const |
void | ComputeTransformation (const G4int copyNo, G4VPhysicalVolume *physVol) const |
void | ComputeDimensions (G4Polyhedra &phedra, const G4int copyNo, const G4VPhysicalVolume *physVol) const |
Definition at line 183 of file G4ParameterisationPolyhedra.hh.
G4ParameterisationPolyhedraZ::G4ParameterisationPolyhedraZ | ( | EAxis | axis, | |
G4int | nCopies, | |||
G4double | offset, | |||
G4double | step, | |||
G4VSolid * | motherSolid, | |||
DivisionType | divType | |||
) |
Definition at line 426 of file G4ParameterisationPolyhedra.cc.
References G4VDivisionParameterisation::CalculateNDiv(), CheckParametersValidity(), DivNDIV, DivWIDTH, G4VDivisionParameterisation::fnDiv, G4VDivisionParameterisation::foffset, G4VDivisionParameterisation::fwidth, G4cout, G4endl, G4PolyhedraHistorical::Num_z_planes, G4VDivisionParameterisation::SetType(), G4VDivisionParameterisation::verbose, and G4PolyhedraHistorical::Z_values.
00429 : G4VParameterisationPolyhedra( axis, nDiv, width, offset, msolid, divType ), 00430 fNSegment(0), 00431 fOrigParamMother(((G4Polyhedra*)fmotherSolid)->GetOriginalParameters()) 00432 { 00433 CheckParametersValidity(); 00434 SetType( "DivisionPolyhedraZ" ); 00435 00436 if( divType == DivWIDTH ) 00437 { 00438 fnDiv = 00439 CalculateNDiv( fOrigParamMother->Z_values[fOrigParamMother->Num_z_planes-1] 00440 - fOrigParamMother->Z_values[0] , width, offset ); 00441 } 00442 else if( divType == DivNDIV ) 00443 { 00444 fwidth = 00445 CalculateNDiv( fOrigParamMother->Z_values[fOrigParamMother->Num_z_planes-1] 00446 - fOrigParamMother->Z_values[0] , nDiv, offset ); 00447 } 00448 00449 #ifdef G4DIVDEBUG 00450 if( verbose >= 1 ) 00451 { 00452 G4cout << " G4ParameterisationPolyhedraZ - # divisions " << fnDiv << " = " 00453 << nDiv << G4endl 00454 << " Offset " << foffset << " = " << offset << G4endl 00455 << " Width " << fwidth << " = " << width << G4endl; 00456 } 00457 #endif 00458 }
G4ParameterisationPolyhedraZ::~G4ParameterisationPolyhedraZ | ( | ) |
void G4ParameterisationPolyhedraZ::CheckParametersValidity | ( | ) | [virtual] |
Reimplemented from G4VDivisionParameterisation.
Definition at line 510 of file G4ParameterisationPolyhedra.cc.
References G4VDivisionParameterisation::CheckParametersValidity(), DivNDIV, DivNDIVandWIDTH, DivWIDTH, FatalException, G4VDivisionParameterisation::fDivisionType, G4VDivisionParameterisation::fmotherSolid, G4VDivisionParameterisation::fnDiv, G4VDivisionParameterisation::foffset, G4VDivisionParameterisation::fReflectedSolid, G4VDivisionParameterisation::fwidth, G4endl, G4Exception(), G4VSolid::GetName(), G4PolyhedraHistorical::Num_z_planes, and G4PolyhedraHistorical::Z_values.
Referenced by G4ParameterisationPolyhedraZ().
00511 { 00512 G4VDivisionParameterisation::CheckParametersValidity(); 00513 00514 // Division will be following the mother polyhedra segments 00515 if( fDivisionType == DivNDIV ) { 00516 if( fOrigParamMother->Num_z_planes-1 != fnDiv ) { 00517 std::ostringstream message; 00518 message << "Configuration not supported." << G4endl 00519 << "Division along Z will be done splitting in the defined" 00520 << G4endl 00521 << "Z planes, i.e, the number of division would be :" 00522 << fOrigParamMother->Num_z_planes-1 << " instead of " 00523 << fnDiv << " !"; 00524 G4Exception("G4ParameterisationPolyhedraZ::CheckParametersValidity()", 00525 "GeomDiv0001", FatalException, message); 00526 } 00527 } 00528 00529 // Division will be done within one polyhedra segment 00530 // with applying given width and offset 00531 if( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH ) { 00532 // Check if divided region does not span over more 00533 // than one z segment 00534 00535 G4int isegstart = -1; // number of the segment containing start position 00536 G4int isegend = -1; // number of the segment containing end position 00537 00538 if ( ! fReflectedSolid ) { 00539 // The start/end position of the divided region 00540 G4double zstart 00541 = fOrigParamMother->Z_values[0] + foffset; 00542 G4double zend 00543 = fOrigParamMother->Z_values[0] + foffset + fnDiv* fwidth; 00544 00545 G4int counter = 0; 00546 while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 ) { 00547 // first segment 00548 if ( zstart >= fOrigParamMother->Z_values[counter] && 00549 zstart < fOrigParamMother->Z_values[counter+1] ) { 00550 isegstart = counter; 00551 } 00552 // last segment 00553 if ( zend > fOrigParamMother->Z_values[counter] && 00554 zend <= fOrigParamMother->Z_values[counter+1] ) { 00555 isegend = counter; 00556 } 00557 ++counter; 00558 } 00559 } 00560 else { 00561 // The start/end position of the divided region 00562 G4double zstart 00563 = fOrigParamMother->Z_values[0] - foffset; 00564 G4double zend 00565 = fOrigParamMother->Z_values[0] - ( foffset + fnDiv* fwidth); 00566 00567 G4int counter = 0; 00568 while ( isegend < 0 && counter < fOrigParamMother->Num_z_planes-1 ) { 00569 // first segment 00570 if ( zstart <= fOrigParamMother->Z_values[counter] && 00571 zstart > fOrigParamMother->Z_values[counter+1] ) { 00572 isegstart = counter; 00573 } 00574 // last segment 00575 if ( zend < fOrigParamMother->Z_values[counter] && 00576 zend >= fOrigParamMother->Z_values[counter+1] ) { 00577 isegend = counter; 00578 } 00579 ++counter; 00580 } 00581 } 00582 00583 if ( isegstart != isegend ) { 00584 std::ostringstream message; 00585 message << "Configuration not supported." << G4endl 00586 << "Division with user defined width." << G4endl 00587 << "Solid " << fmotherSolid->GetName() << G4endl 00588 << "Divided region is not between two Z planes."; 00589 G4Exception("G4ParameterisationPolyhedraZ::CheckParametersValidity()", 00590 "GeomDiv0001", FatalException, message); 00591 } 00592 00593 fNSegment = isegstart; 00594 } 00595 }
void G4ParameterisationPolyhedraZ::ComputeDimensions | ( | G4Polyhedra & | phedra, | |
const G4int | copyNo, | |||
const G4VPhysicalVolume * | physVol | |||
) | const [virtual] |
Reimplemented from G4VPVParameterisation.
Definition at line 649 of file G4ParameterisationPolyhedra.cc.
References DivNDIV, DivNDIVandWIDTH, DivWIDTH, G4VSolid::DumpInfo(), G4VDivisionParameterisation::fDivisionType, G4VDivisionParameterisation::foffset, G4VDivisionParameterisation::fReflectedSolid, G4VDivisionParameterisation::fwidth, G4cout, G4endl, G4PolyhedraHistorical::Num_z_planes, G4PolyhedraHistorical::numSide, G4PolyhedraHistorical::Opening_angle, G4Polyhedra::Reset(), G4PolyhedraHistorical::Rmax, G4PolyhedraHistorical::Rmin, G4Polyhedra::SetOriginalParameters(), G4PolyhedraHistorical::Start_angle, G4VDivisionParameterisation::verbose, and G4PolyhedraHistorical::Z_values.
00651 { 00652 // Define division solid 00653 G4PolyhedraHistorical origparam; 00654 G4int nz = 2; 00655 origparam.Num_z_planes = nz; 00656 origparam.numSide = fOrigParamMother->numSide; 00657 origparam.Start_angle = fOrigParamMother->Start_angle; 00658 origparam.Opening_angle = fOrigParamMother->Opening_angle; 00659 00660 // Define division solid z sections 00661 origparam.Z_values = new G4double[nz]; 00662 origparam.Rmin = new G4double[nz]; 00663 origparam.Rmax = new G4double[nz]; 00664 origparam.Z_values[0] = - fwidth/2.; 00665 origparam.Z_values[1] = fwidth/2.; 00666 00667 if ( fDivisionType == DivNDIV ) { 00668 // The position of the centre of copyNo-th mother polycone segment 00669 G4double posi = ( fOrigParamMother->Z_values[copyNo] 00670 + fOrigParamMother->Z_values[copyNo+1])/2; 00671 00672 origparam.Z_values[0] = fOrigParamMother->Z_values[copyNo] - posi; 00673 origparam.Z_values[1] = fOrigParamMother->Z_values[copyNo+1] - posi; 00674 origparam.Rmin[0] = fOrigParamMother->Rmin[copyNo]; 00675 origparam.Rmin[1] = fOrigParamMother->Rmin[copyNo+1]; 00676 origparam.Rmax[0] = fOrigParamMother->Rmax[copyNo]; 00677 origparam.Rmax[1] = fOrigParamMother->Rmax[copyNo+1]; 00678 } 00679 00680 if ( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH ) { 00681 if ( ! fReflectedSolid ) { 00682 origparam.Z_values[0] = - fwidth/2.; 00683 origparam.Z_values[1] = fwidth/2.; 00684 00685 // The position of the centre of copyNo-th division 00686 G4double posi 00687 = fOrigParamMother->Z_values[0] + foffset + (2*copyNo + 1) * fwidth/2.; 00688 00689 // The first and last z sides z values 00690 G4double zstart = posi - fwidth/2.; 00691 G4double zend = posi + fwidth/2.; 00692 origparam.Rmin[0] = GetRmin(zstart, fNSegment); 00693 origparam.Rmax[0] = GetRmax(zstart, fNSegment); 00694 origparam.Rmin[1] = GetRmin(zend, fNSegment); 00695 origparam.Rmax[1] = GetRmax(zend, fNSegment); 00696 } 00697 else { 00698 origparam.Z_values[0] = fwidth/2.; 00699 origparam.Z_values[1] = - fwidth/2.; 00700 00701 // The position of the centre of copyNo-th division 00702 G4double posi 00703 = fOrigParamMother->Z_values[0] - ( foffset + (2*copyNo + 1) * fwidth/2.); 00704 00705 // The first and last z sides z values 00706 G4double zstart = posi + fwidth/2.; 00707 G4double zend = posi - fwidth/2.; 00708 origparam.Rmin[0] = GetRmin(zstart, fNSegment); 00709 origparam.Rmax[0] = GetRmax(zstart, fNSegment); 00710 origparam.Rmin[1] = GetRmin(zend, fNSegment); 00711 origparam.Rmax[1] = GetRmax(zend, fNSegment); 00712 } 00713 00714 // It can happen due to rounding errors 00715 if ( origparam.Rmin[0] < 0.0 ) origparam.Rmin[0] = 0.0; 00716 if ( origparam.Rmin[nz-1] < 0.0 ) origparam.Rmin[1] = 0.0; 00717 } 00718 00719 phedra.SetOriginalParameters(&origparam); // copy values & transfer pointers 00720 phedra.Reset(); // reset to new solid parameters 00721 00722 #ifdef G4DIVDEBUG 00723 if( verbose >= 2 ) 00724 { 00725 G4cout << "G4ParameterisationPolyhedraZ::ComputeDimensions()" << G4endl 00726 << "-- Parametrised phedra copy-number: " << copyNo << G4endl; 00727 phedra.DumpInfo(); 00728 } 00729 #endif 00730 }
void G4ParameterisationPolyhedraZ::ComputeTransformation | ( | const G4int | copyNo, | |
G4VPhysicalVolume * | physVol | |||
) | const [virtual] |
Implements G4VDivisionParameterisation.
Definition at line 600 of file G4ParameterisationPolyhedra.cc.
References G4VDivisionParameterisation::ChangeRotMatrix(), DivNDIV, DivNDIVandWIDTH, DivWIDTH, G4VDivisionParameterisation::faxis, G4VDivisionParameterisation::fDivisionType, G4VDivisionParameterisation::foffset, G4VDivisionParameterisation::fReflectedSolid, G4VDivisionParameterisation::fwidth, G4cout, G4endl, G4VPhysicalVolume::SetTranslation(), G4VDivisionParameterisation::verbose, and G4PolyhedraHistorical::Z_values.
00601 { 00602 if ( fDivisionType == DivNDIV ) { 00603 // The position of the centre of copyNo-th mother polycone segment 00604 G4double posi = ( fOrigParamMother->Z_values[copyNo] 00605 + fOrigParamMother->Z_values[copyNo+1])/2; 00606 physVol->SetTranslation( G4ThreeVector(0, 0, posi) ); 00607 } 00608 00609 if ( fDivisionType == DivNDIVandWIDTH || fDivisionType == DivWIDTH ) { 00610 // The position of the centre of copyNo-th division 00611 00612 G4double posi = fOrigParamMother->Z_values[0]; 00613 00614 if ( ! fReflectedSolid ) 00615 posi += foffset + (2*copyNo + 1) * fwidth/2.; 00616 else 00617 posi -= foffset + (2*copyNo + 1) * fwidth/2.; 00618 00619 physVol->SetTranslation( G4ThreeVector(0, 0, posi) ); 00620 } 00621 00622 //----- calculate rotation matrix: unit 00623 00624 #ifdef G4DIVDEBUG 00625 if( verbose >= 2 ) 00626 { 00627 G4cout << " G4ParameterisationPolyhedraZ - position: " << posi << G4endl 00628 << " copyNo: " << copyNo << " - foffset: " << foffset/deg 00629 << " - fwidth: " << fwidth/deg << G4endl; 00630 } 00631 #endif 00632 00633 ChangeRotMatrix( physVol ); 00634 00635 #ifdef G4DIVDEBUG 00636 if( verbose >= 2 ) 00637 { 00638 G4cout << std::setprecision(8) << " G4ParameterisationPolyhedraZ " 00639 << copyNo << G4endl 00640 << " Position: " << origin << " - Width: " << fwidth 00641 << " - Axis: " << faxis << G4endl; 00642 } 00643 #endif 00644 }
G4double G4ParameterisationPolyhedraZ::GetMaxParameter | ( | ) | const [virtual] |
Implements G4VDivisionParameterisation.
Definition at line 503 of file G4ParameterisationPolyhedra.cc.
References G4PolyhedraHistorical::Num_z_planes, and G4PolyhedraHistorical::Z_values.
00504 { 00505 return std::abs (fOrigParamMother->Z_values[fOrigParamMother->Num_z_planes-1] 00506 -fOrigParamMother->Z_values[0]); 00507 }