spradlin / WCSim

The WCSim GEANT4 application
0 stars 0 forks source link

Hunting overlaps #24

Closed spradlin closed 1 year ago

spradlin commented 1 year ago

Using the developments described in #20 and #23, search for (and resolve) overlaps in the HyperKWithOD detector geometry.

spradlin commented 1 year ago

Failed overlap checks

Looking for failed overlap checks in the construction phase. I am using the same test scripts as i did in #15. The only change that i have made is to reduce the number of generated events to 1. I ran with the overlap checking turned on---it produces a 233297-line log.

Successful overlap checks produce log messages of the form

Checking overlaps for volume <element name> ... OK!

So i did a preliminary search for failed overlap checks that do not have 'OK!' with the following command:

$ grep "Checking overlaps for volume" wcsim_hkod_explicit_SensitiveDetector_Only_NoTrigger_gun5000MeVmumOD_dark0_n1_seed98.overlap.log | grep -v OK
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCPMT ... 
Checking overlaps for volume WCPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCMultiPMT ... 
Checking overlaps for volume WCPMT ... 
Checking overlaps for volume WCPMT ... 
Checking overlaps for volume TopCapAssembly ... 
Checking overlaps for volume BottomCapAssembly ... 

So, 14 potential failures of the overlap check on placement. Some of these might be in the initial SuperK construction. Let me try to track these down in the log file and get some context.

spradlin commented 1 year ago

Overlaps detected during SuperK construction

Just clip some context around the failed overlap checks.

1: SuperK WCMultiPMT

total on cap: 1672
Coverage was calculated to be: 0.380348
Debug B.Q : place border pmt, spacing horiz = 707.079, vertical = -699.843
Add PMT on barrel cell 0, 0
Checking overlaps for volume WCMultiPMT ... 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16727.6,-984.589,942.197), overlapping by at least: 3.94316 mm 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

2: Super WCMultiPMT

Add PMT on barrel cell 1, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16660.8,-331.933,937.869), overlapping by at least: 4.81568 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

3: SuperK WCMultiPMT

Add PMT on barrel cell 2, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16677.2,260.669,936.256), overlapping by at least: 3.53979 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

4: SuperK WCMultiPMT

Add PMT on barrel cell 3, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16660.8,984.5,877.626), overlapping by at least: 5.55886 mm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

5: SuperK WCPMT

Debug B.Q : Add extra tower cap
Add PMTs in extra tower, cell 0, 0
Checking overlaps for volume WCPMT ... 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCPMT
          with its mother volume WCspecialBarrelBorderCell
          at mother local point (16716.5,-335.171,953.578), overlapping by at least: 5.14481 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

6: SuperK WCPMT

Add PMTs in extra tower, cell 1, 0
Checking overlaps for volume WCPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCPMT
          with its mother volume WCspecialBarrelBorderCell
          at mother local point (16671,-1164.76,921.1), overlapping by at least: 3.60122 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

7: SuperK WCMultiPMT

Debug B.Q : place border pmt, spacing horiz = 707.079, vertical = 699.843
Add PMT on barrel cell 0, 0
Checking overlaps for volume WCMultiPMT ... 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16662,-997.605,-945.903), overlapping by at least: 5.2945 cm 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

8: SuperK WCMultiPMT

Add PMT on barrel cell 1, 0
Checking overlaps for volume WCMultiPMT ... 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16699.1,-454.423,-932.949), overlapping by at least: 1.75268 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

9: SuperK WCMultiPMT

Add PMT on barrel cell 2, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16712.5,374.126,-953.007), overlapping by at least: 2.22681 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

10: SuperK WCMultiPMT

Add PMT on barrel cell 3, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (16731.5,1124.47,-945.687), overlapping by at least: 3.6801 mm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

11: SuperK WCPMT

Debug B.Q : Add extra tower cap
Add PMTs in extra tower, cell 0, 0
Checking overlaps for volume WCPMT ... 
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCPMT
          with its mother volume WCspecialBarrelBorderCell
          at mother local point (16741.4,-224.563,-921.99), overlapping by at least: 1.47832 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

12: SuperK WCPMT

Add PMTs in extra tower, cell 1, 0
Checking overlaps for volume WCPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCPMT
          with its mother volume WCspecialBarrelBorderCell
          at mother local point (16697.1,-921.345,-916.725), overlapping by at least: 2.1716 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

13: SuperK TopCapAssembly

Checking overlaps for volume TopCapAssembly ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume TopCapAssembly
          with its mother volume WCBarrel
          at mother local point (-5549.3,9201.77,18121), overlapping by at least: 1 mm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

14: SuperK BottomCapAssembly

Checking overlaps for volume BottomCapAssembly ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume BottomCapAssembly
          with its mother volume WCBarrel
          at mother local point (3427.36,-9267.83,-18121), overlapping by at least: 1 mm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------
spradlin commented 1 year ago

Recheck SuperK

From the previous, it looks like all of the detected overlaps occur happen in the default SuperK construction phase. Wait, is SuperK still the default? Yes, according to the constructor WCSimDetectorConstruction::WCSimDetectorConstruction().

OK, check for overlaps in the default SuperK by running a macro that does not perform a second detector construction and just uses the default. The steering macro for this test has just a single line:

/run/beamOn 1

All 14 of the previously recorded overlaps are present in this test. All of the overlaps appear to occur in the SuperK geometry.

OK, well, that's a place to start.

spradlin commented 1 year ago

The WCBarrelBorderCell

Many of the overlap errors are reported for placement of PMTs or MultiPMTs in WCBarrelBorderCell, and specifically for the first column of PMTs in those cells (vertical index 0). This looks like the horizontal offset of the first column may be wrong.

The PMTs are placed in the logicWCBarrelBorderCell in WCSimDetectorConstruction::ConstructCaps() L2010-L2021. It places one of two logical PMT composite volumes, depending on a string of conditions:

(i==j && hybrid && WCPMTPercentCoverage2!=0)?logicWCPMT2:logicWCPMT,

There are a few things problematic with this:

Anyway, the conditional should evaluate to false for the initial SuperK geometry becuase of the hybrid default. This means the second logical PMT, logicWCPMT should be the one that is placed.

logicWCPMT is a local variable declared and initialized in L1875-L1877:

  G4LogicalVolume* logicWCPMT;
  if(nID_PMTs<=1) logicWCPMT = ConstructPMT(WCPMTName, WCIDCollectionName,"tank",nID_PMTs);
  else logicWCPMT = ConstructMultiPMT(WCPMTName, WCIDCollectionName,"tank",nID_PMTs);

The SuperK geometry calls WCSimDetectorConstruction::InitSinglePMT() which sets WCSimDetectorConstruction::nID_PMTs to 1. So logicWCPMT is constructed by one of the two overloaded instances of WCSimDetectorConstruction::ConstructPMT(). Since WCSimDetectorConstruction::nID_PMTs is a G4int, it should be the one defined in WCSimConstructPMT.cc L32-L478 .

The fact that this snippet calls WCSimDetectorConstruction::ConstructPMT() if nID_PMTs==1 and WCSimDetectorConstruction::ConstructMultiPMT() otherwise prompts the question of whether WCSimDetectorConstruction::ConstructPMT() is ever called with nID_PMTs>1. It has some code that is invoked in this case, but does any client call it this way? Search for invocations of it:

Ok, these last two are instances where the MultiPMT builder calls the single PMT builder.

spradlin commented 1 year ago

All of the visualization configuration might be better served by a configurable factory than by these hard-coded things.

spradlin commented 1 year ago

Factorizing the PMT construction

I believe that many of the problems in the geometry construction come from its monolithic construction with a lot of convoluted options. This makes consistent dimension checking nearly impossible.

One strategy that might help is further factorization of the construction into components, if it is implemented well. The components can be represented by classes or objects that can be queried for dimension information. This kind of encapsulation also usually makes debugging easier.

The PMT logical volume consists of

One layer of recursion is allowed for multiPMTs. The containing volume is the most important for interfaces with the rest of the geometry.

spradlin commented 1 year ago

PMT containing volumes

Single PMT version 1

This is the version created by the nuPRISM version of the overloaded WCSimDetectorConstruction::ConstructPMT(). In the current code, it seems to be used for all single PMTs that are not accompanied by a wavelength-shifting plate.

Its solid implementation is a simple solid cylinder, though implemented as a Polycon (not sure why). There are two variants:

  1. No reflector---dimensions are determined by only the WCSimPMTObject data
  2. With reflector---dimensions are determined by a combination of the WCSimPMTObject data and the dimensions of the reflector.

Single PMT within a multiPMT

Also created in the nuPRISM version of the overloaded WCSimDetectorConstruction::ConstructPMT().

Its solid implementation is a segment of a spherical shell.

Come back to this.

spradlin commented 1 year ago

Visualization of the SuperK WCBarrelBorderCell

SuperK_WCBarrelBorderCell_0000

I was able to make this visualization of the SuperK WCBarrelBorderCell. It looks like the cylindrical PMT volumes along the beveled edge jut out of the bevel. Probably need to increase the margin on that edge in order to prevent that.

spradlin commented 1 year ago

WCBarrelBorderCell construction

There is probably a straightforward way to fix the margin, but i need to know a little more about the WCBarrelBorderCell volume.

Its solid implementation is a G4Polyhedra:

  1496    G4Polyhedra* solidWCBarrelBorderCell = new G4Polyhedra("WCBarrelBorderCell",  
  1497                                                                                                                   -dPhi/2., // phi start
  1498                                                                                                                   dPhi, //total angle
  1499                                                                                                                   1, //NPhi-gon
  1500                                                                                                                   3,
  1501                                                                                                                   borderAnnulusZ,
  1502                                                                                                                   borderAnnulusRmin,
  1503                                                                                                                   borderAnnulusRmax);

The arrays that define polyhedra is in L1447-L1457:

    39  #define MIRROR_WCSIM_DEVELOP_POLY
     :
  1447    G4double borderAnnulusZ[3] = {
  1448  #ifdef MIRROR_WCSIM_DEVELOP_POLY
  1449          -barrelCellHeight/2.*zflip,
  1450          (-barrelCellHeight/2.+(WCIDRadius-innerAnnulusRadius))*zflip,
  1451  #else
  1452          (-barrelCellHeight/2.-(WCIDRadius-innerAnnulusRadius))*zflip,
  1453          -barrelCellHeight/2.*zflip,
  1454  #endif
  1455          barrelCellHeight/2.*zflip};
  1456    G4double borderAnnulusRmin[3] = { WCIDRadius, innerAnnulusRadius, innerAnnulusRadius};
  1457    G4double borderAnnulusRmax[3] = {outerAnnulusRadius, outerAnnulusRadius,outerAnnulusRadius};

(The conditional compilation can't be good.)

I think that the bevel is between the first two z-positions. With that assumption, it looks like there should be an extra WCIDRadius-innerAnnulusRadius margin.

OK, let us try to track these dimensions. The fundamental dimensions defined in the detector config are

The WCIDRadius is computed from WCIDDiameter

    69  G4LogicalVolume* WCSimDetectorConstruction::ConstructCylinder()
    70  {
     :
    81    WCIDRadius = WCIDDiameter/2.;
    :
  1420  }

The innerAnnulusRadius has a conditional value

    69  G4LogicalVolume* WCSimDetectorConstruction::ConstructCylinder()
    70  {
     :
   100    if(hybrid){
   101          if(mPMT_vessel_cyl_height + mPMT_vessel_radius < 1.*mm)
   102            innerAnnulusRadius = WCIDRadius - std::max(WCPMTExposeHeight,WCPMTExposeHeight2)-1.*mm;
   103          else
   104            innerAnnulusRadius = WCIDRadius - std::max(WCPMTExposeHeight,(mPMT_vessel_cyl_height + mPMT_vessel_radius)) -1.*mm;
   105    }
   106    else{
   107          if(mPMT_vessel_cyl_height + mPMT_vessel_radius < 1.*mm)
   108            innerAnnulusRadius = WCIDRadius - WCPMTExposeHeight -1.*mm;
   109          else
   110            innerAnnulusRadius = WCIDRadius - (mPMT_vessel_cyl_height + mPMT_vessel_radius) -1.*mm;
   111    }
    :
  1420  }

The set of values defined for the SuperK geometry put us in the assignment defined on L108. The innerAnnulusRadius is smaller than the WCIDRadius by a bit more than the expose height of the PMT.

Because the SuperK does not define an OD, the definion of outerAnnulusRadius on L114 applies

    69  G4LogicalVolume* WCSimDetectorConstruction::ConstructCylinder()
    70  {
     :
   114    outerAnnulusRadius = WCIDRadius + WCBlackSheetThickness + 1.*mm;//+ Stealstructure etc.
   115    if(isODConstructed){
   116      const G4double sphereRadius =
   117            (WCPMTODExposeHeight*WCPMTODExposeHeight+ WCPMTODRadius*WCPMTODRadius)/(2*WCPMTODExposeHeight);
   118      WCODRadius = 
   119            WCIDRadius + WCBlackSheetThickness + WCODDeadSpace + // ID Structure  
   120            WCODTyvekSheetThickness;  // Tyvek attached to structure
   121      outerAnnulusRadius = WCODRadius + sphereRadius;                        122    }
    :
  1420  }

So the following heirarchy applies:

innerAnnulusRadius < WCIDRadius < outerAnnulusRadius

The thickness of the bulk of the WCBarrelBorderCell is

outerAnnulusRadius - innerAnnulusRadius = WCPMTExposeHeight + WCBlackSheetThickness + 2mm

which the bevel narrows to

outerAnnulusRadius - WCIDRadius = WCBlackSheetThickness + 1mm

The width of the bevel in $z$ is the difference between the first two elements of borderAnnulusZ

bevel = WCIDRadius - innerAnnulusRadius = WCPMTExposeHeight + 1mm

spradlin commented 1 year ago

PMT placement in WCBarrelBorderCell

Let's look at the loop that places PMTs in the WCBarrelBorderCell

  1423  G4LogicalVolume* WCSimDetectorConstruction::ConstructCaps(G4int zflip)
  1424  {
     :
  2004           
  2005          for(G4double i = 0; i < WCPMTperCellHorizontal; i++){
  2006            for(G4double j = 0; j < WCPMTperCellVertical; j++){
  2007                  G4ThreeVector PMTPosition =  G4ThreeVector(WCIDRadius,
  2008                                                                                                     -barrelCellWidth/2.+(i+0.5)*horizontalSpacing,
  2009                                                                                                     (-(barrelCellHeight-WCBorderPMTOffset)/2.+(j+0.5)*verticalSpacing)*zflip);     
  2010  #ifdef ACTIVATE_IDPMTS                            
  2011  #ifdef WCSIMCONSTRUCTCYLINDER_VERBOSE                   
  2012                  G4cout << "Add PMT on barrel cell " << i << ", " << j << G4endl;
  2013  #endif
  2014                  //G4VPhysicalVolume* physiWCBarrelBorderPMT =
  2015                  new G4PVPlacement(WCPMTRotation,                      // its rotation                                                   
  2016                                                    PMTPosition,
  2017                                                    (i==j && hybrid && WCPMTPercentCoverage2!=0)?logicWCPMT2:logicWCPMT,                // its logical volume
  2018                                                    pmtname,             // its name
  2019                                                    logicWCBarrelBorderCell,         // its mother volume
  2020                                                    false,                     // no boolean operations
  2021                                                    (int)(i*WCPMTperCellVertical+j),
  2022                                                    checkOverlapsPMT);
  2023  #endif
  2024
  2025                  // logicWCPMT->GetDaughter(0),physiCapPMT is the glass face. If you add more
  2026                  // daugter volumes to the PMTs (e.g. a acryl cover) you have to check, if
  2027                  // this is still the case.
  2028            }
  2029          }
     :
  2199  }

There is this WCBorderPMTOffset in the $z$ coordinate that might be there to prevent this kind of overlap. Maybe something is going wrong with it. Where does it come from? It is initialized to 0 in the WCSimDetectorConstruction constructor and then reassigned in two geometry configurations---HyperKWithOD and HyperK_HybridmPMT_WithOD---but not in SuperK.

OK, so it looks like a more-or-less correct solution is to update the value of WCBorderPMTOffset when the dimensions of the WCBarrelBorderRing are defined. I also needed to fix the arithmetic in the placement position---the offset should not be divided by 2 in the placement calculations.

That seems to fix the WCPMT placement overlaps. Let me get a clean branch with the fix for a pull request. The endcap assembly overlaps with the barrel remain.

spradlin commented 1 year ago

SuperK endcap overlaps

The remaining overlaps identified by GEANT4 overlap checking occur when the endcap assemblies are placed in the WCBarrel:

Checking overlaps for volume TopCapAssembly ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume TopCapAssembly
          with its mother volume WCBarrel
          at mother local point (11387.8,-2038.54,18121), overlapping by at least: 1 mm 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

Checking overlaps for volume BottomCapAssembly ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume BottomCapAssembly
          with its mother volume WCBarrel
          at mother local point (8571.87,14312.8,-18121), overlapping by at least: 1 mm 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

The size of the overlap, 1 mm, is highly suggestive. There are several places in which 1 mm margins are introduced, probably to avoid overlaps.

spradlin commented 1 year ago

WCBarrel construction

The WCBarrel is a cylinder:

    69  G4LogicalVolume* WCSimDetectorConstruction::ConstructCylinder()
    70  {
     :
   216    //-----------------------------------------------------
   217    // everything else is contained in this water tubs
   218    //-----------------------------------------------------
   219    G4Tubs* solidWCBarrel = new G4Tubs("WCBarrel",
   220                                                                           0.0*m, 
   221                                                                           WCRadius+1.*m, // add a bit of extra space
   222                                                                           .5*WCLength,  //jl145 - per blueprint
   223                                                                           0.*deg,
   224                                                                           360.*deg);
     :    
   228    G4LogicalVolume* logicWCBarrel =              
   229      new G4LogicalVolume(solidWCBarrel,          
   230                                                  G4Material::GetMaterial(water), 
   231                                                  "WCBarrel",
   232                                                  0,0,0);
     :
  1420  }

Fundamental dimensions from the config

The dimensions of WCBarrel are computed in ConstructCylinder():

    69  G4LogicalVolume* WCSimDetectorConstruction::ConstructCylinder()
    70  {
     :  
    81    WCIDRadius = WCIDDiameter/2.;
     :
    83    WCBarrelRingNPhi     = (G4int)(WCBarrelNumPMTHorizontal/WCPMTperCellHorizontal); 
     :
    85    totalAngle  = 2.0*pi*rad*(WCBarrelRingNPhi*WCPMTperCellHorizontal/WCBarrelNumPMTHorizontal) ;
     :
    87    dPhi        =  totalAngle/ WCBarrelRingNPhi;
     :
    89    barrelCellHeight  = (WCIDHeight-2.*WCBarrelPMTOffset)/WCBarrelNRings;
     :
    91    mainAnnulusHeight = WCIDHeight -2.*WCBarrelPMTOffset -2.*barrelCellHeight;
     :
   113    //TF: need to add a Polyhedra on the other side of the outerAnnulusRadius for the OD
   114    outerAnnulusRadius = WCIDRadius + WCBlackSheetThickness + 1.*mm;//+ Stealstructure etc.
   115    if(isODConstructed){
   116      const G4double sphereRadius =
   117            (WCPMTODExposeHeight*WCPMTODExposeHeight+ WCPMTODRadius*WCPMTODRadius)/(2*WCPMTODExposeHeight);
   118      WCODRadius =
   119            WCIDRadius + WCBlackSheetThickness + WCODDeadSpace + // ID Structure
   120            WCODTyvekSheetThickness;  // Tyvek attached to structure
   121      outerAnnulusRadius = WCODRadius + sphereRadius;
   122    }
     :
   147    WCLength    = WCIDHeight + 2*(WCODHeightWaterDepth + WCBlackSheetThickness + WCODDeadSpace + WCODTyvekSheetThickness);
   148    WCRadius    = (outerAnnulusRadius + WCODLateralWaterDepth)/cos(dPhi/2.) ; 
   149    if(isODConstructed){
   150      WCLength    = WCIDHeight + 2*(WCODHeightWaterDepth + WCBlackSheetThickness + WCODDeadSpace + WCODTyvekSheetThickness);
   151      WCRadius    = (outerAnnulusRadius + WCODLateralWaterDepth)/cos(dPhi/2.) ;
   152    }
     :
  1420  }

Note: The conditional clause that begins on L149 is redundant. There is no difference in the calculation with and without isODConstructed. I am not sure that this is intentional.

Endcap construction

The endcap assembly volumes are constructed as cylinders:

  1423  G4LogicalVolume* WCSimDetectorConstruction::ConstructCaps(G4int zflip)
  1424  {
  1425    capAssemblyHeight = (WCIDHeight-mainAnnulusHeight)/2+1*mm+WCBlackSheetThickness;
  1426  
  1427    G4Tubs* solidCapAssembly = new G4Tubs("CapAssembly",
  1428                                                                                  0.0*m,
  1429                                                                                  outerAnnulusRadius/cos(dPhi/2.), 
  1430                                                                                  capAssemblyHeight/2,
  1431                                                                                  0.*deg,
  1432                                                                                  360.*deg);
  1433
  1434    G4LogicalVolume* logicCapAssembly =
  1435      new G4LogicalVolume(solidCapAssembly,
  1436                          G4Material::GetMaterial(water),
  1437                          "CapAssembly",
  1438                          0,0,0);
     : 
  2197    return logicCapAssembly;
  2198 
  2199  }

The difference between the radius of thecontaining volume and the radius of the endcap assembly is

WCRadius - outerAnnulusRadius/cos(dPhi/2.) = WCODLateralWaterDepth/cos(dPhi/2.)

As long as WCODLateralWaterDepth >= 0, we would expect the endcap assembly to fit radially. Unfortunately, WCODLateralWaterDepth is only initialized in the HyperKWithOD and HyperK_HybridmPMT_WithOD conigurations. It is possible to set it with a macro command and it defaults to 1 (i am not certain about the units). However, if the command is not called, it is UNINITIALIZED.

So, the radius should fit.

Endcap placement

The endcap assemblies are the last volumes placed in the WCBarrel:

    69  G4LogicalVolume* WCSimDetectorConstruction::ConstructCylinder()
    70  {  
     :          
  1401    //G4VPhysicalVolume* physiTopCapAssembly =
  1402          new G4PVPlacement(0,
  1403                                            G4ThreeVector(0.,0.,(mainAnnulusHeight/2.+ capAssemblyHeight/2.)),
  1404                                            logicTopCapAssembly,
  1405                                            "TopCapAssembly",
  1406                                            logicWCBarrel,
  1407                                            false, 0,
  1408                                            checkOverlaps);
  1409 
  1410    //G4VPhysicalVolume* physiBottomCapAssembly =
  1411      new G4PVPlacement(0,
  1412                                            G4ThreeVector(0.,0.,(-mainAnnulusHeight/2.- capAssemblyHeight/2.)),
  1413                                            logicBottomCapAssembly,
  1414                                            "BottomCapAssembly",
  1415                                            logicWCBarrel,
  1416                                            false, 0,
  1417                                            checkOverlaps);
  1418 
  1419    return logicWC;
  1420  }

So, the requirement for the endcap assembly to fit longitudinally in the barrel is

WCLength/2 >= mainAnnulusHeight/2 +capAssemblyHeight
(WCIDHeight/2 + WCODHeightWaterDepth + WCBlackSheetThickness + WCODDeadSpace + WCODTyvekSheetThickness) >= mainAnnulusHeight/2 + ((WCIDHeight - mainAnnulusHeight)/2 + 1mm + WCBlackSheetThickness)
>= WCIDHeight/2 + 1
mm + WCBlackSheetThickness
WCODHeightWaterDepth + WCODDeadSpace + WCODTyvekSheetThickness >= 1*mm

So, let us track down what value these quantities have when there is no OD defined. It looks like these are in the same situatation as WCODHeightWaterDepth---uninitialized in most cases. The verbose reporting from a test job indicates that these have value 0:

HYBRID = 0
GEOMCHECK3 WCIDRadius   16840.8
GEOMCHECK3 WCBlackSheetThickness        20
GEOMCHECK3 WCODDeadSpace        0
GEOMCHECK3 WCODTyvekSheetThickness      0
GEOMCHECK3 sphereRadius         16861.8
GEOMCHECK3 WCPMTODExposeHeight  0
GEOMCHECK3 WCPMTODExposeHeight  0
GEOMCHECK3 WCPMTODRadius        0
GEOMCHECK3 WCPMTODRadius        0
GEOMCHECK3 WCPMTODExposeHeight  0
GEOMCHECK2 WCLength     36240 
GEOMCHECK2 WCIDHeight   36200
GEOMCHECK2 WCODHeightWaterDepth         0
GEOMCHECK2 WCBlackSheetThickness        20
GEOMCHECK2 WCODDeadSpace        0
GEOMCHECK2 WCODTyvekSheetThickness      0
GEOMCHECK2 WCRadius     16921.1
GEOMCHECK2 outerAnnulusRadius   16861.8
GEOMCHECK2 WCODLateralWaterDepth        0

Right, so, in order to avoid the endcap assembly sticking out beyond the WCBarrel containing volume, the 1mm fudge factor for capAssemblyHeight must also be added to the WCLength for each endcap.

spradlin commented 1 year ago

Other geometries to check

According to Tom Dealtry's WCSim status report at the HK Collab meeting, HK is validating with the following geometries:

So, let us run overlap checks on these geometries, now that the SuperK geometry is not showing any obvious overlaps.

HyperK20pcWithOD is on the new list of obsolete geometries. Put off checking that geometry.

spradlin commented 1 year ago

HyperK_HybridmPMT overlap check

Starting from my local branch with the fixes identified above for the SuperK overlaps, run 1 event with overlap checking for the HyperK_HybridmPMT geometry.

The following two overlaps were identified:

0, -40, 40, 0, 40, 0, 304
Min viewing angle eta from input file for mPMT filling is : 0.174533
Checking overlaps for parameterised volume pmt ...  
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVParameterised::CheckOverlaps()
Overlap with mother volume !
         Overlap is detected for volume pmt, parameterised instance: 0
          with its mother volume WCMultiPMT_container
          at mother local point (127.754,-0.0324057,33.3208), overlapping by at least: 3.16619 mm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

Checking overlaps for parameterised volume pmt ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVParameterised::CheckOverlaps()
Overlap with mother volume ! 
         Overlap is detected for volume pmt, parameterised instance: 0
          with its mother volume WCMultiPMT_container
          at mother local point (149.809,-16.3052,16.882), overlapping by at least: 5.62 mm 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

Hemisphere has overlapping PMTs: Retry

These overlaps look like they occur in the assembly of a multiPMT.

The last line is interesting. What does it mean? Here is the code fragment that produces this message:

    37  G4LogicalVolume* WCSimDetectorConstruction::ConstructMultiPMT(G4String PMTName, G4String CollectionName, G4String detectorElement, G4int nIDPMTs)
    38  {
     :
   525    // Need a 4mm tolerance : still not perfect though because overlaps with mother volume happen still..
   526    if(hemisphere->CheckOverlaps(1000,4,1)){
   527      G4cerr << "Hemisphere has overlapping PMTs: Retry" << G4endl;
   528      //return NULL;
   529    }
     :
   692  }

It used to be part of a fail state. Now it does nothing but report. The conditional does not appear anywhere else in the code. OK, so nothing is retried.

spradlin commented 1 year ago

Some visulazations of the multiPMT

WCPMT_all_0000

The overall containing volume of the multiPMT is a cylinder.

Because the reported overlaps occur between the volume named WCMultiPMT_container and two of the 'pmt' volumes that are placed in it, let us look specifically at those.

WCPMT_container_all_0001

The yellow volume is WCMultiPMT_container, the blue volumes are the 'pmt's.

WCPMT_container_zoom_0000

A closer look at the interior of WCMultiPMT_container appears to show the blue pmt plugs extending beyond the yellow shell volume that should contain them.

spradlin commented 1 year ago

PMT placement in multiPMT

This occurs in WCSimDetectorConstruction::ConstructMultiPMT(), naturally. The PMT logical volume is referenced by the pointer logicWCPMT; the WCMultiPMT_container into which they are placed is referenced by the pointer logic_mPMT_container

PMT construction

The PMT geometry construction is performed by a [call]() to the non-bool version of WCSimDetectorConstruction::ConstructPMT(https://github.com/tdealtry/WCSim/blob/aa394f6a77956564f8965b4784c3b1944c5f534e/src/WCSimConstructMultiPMT.cc#L417). The last argument indicates which of the two PMT containing volumes to use---the multiPMT version, in the present case.

It is a segment of a spherical shell:

   146      // for single PMTs no overlap issues and will rest on flat black sheets
   147      // Need spherical shell for curved top to better match (and reduce overlaps) with
   148      // gel layer ("container") and mPMT_vessel.
   149      solidWCPMT =
   150        new G4Sphere("WCPMT",
   151                     mPMT_vessel_radius_curv - mPMT_outer_material_d - pmtModuleHeight, //rMin = 342 - 10 - 59.62 = 272.38 mm, 59.62mm is position of support structure
   152                     mPMT_vessel_radius_curv - mPMT_outer_material_d,          //rMax = 332 mm
   153                     0.0*deg,                                             //phiStart
   154                     360.0*deg,                                           //Deltaphi
   155                     0.0*deg,                                             //thetaStart
   156                     mPMT_pmt_openingAngle);                                            //Deltatheta 8.3 deg                                            //Deltatheta

The dimensions WCSimDetectorConstruction::mPMT_vessel_radius_curv and WCSimDetectorConstruction::mPMT_outer_material_d are defined in the detector configuration. For the HyperK_HybridmPMT geometry,

   528  void WCSimDetectorConstruction::SetHyperK_HybridmPMTGeometry()
   529  { 
     :
   550    //mPMT params go first because detector depends on it:
   551    mPMT_vessel_cyl_height = 38.*CLHEP::mm;    //option A, option B would be 277 mm
   552    mPMT_vessel_radius_curv = 342.*CLHEP::mm;  //needs to include the vessel thickness, as we construct from outside inwards.
   553    mPMT_vessel_radius = 254.*CLHEP::mm;
     :
   560    mPMT_outer_material_d = 10*CLHEP::mm;
     :
   592  } 

The value of pmtModuleHeight is a HARD-CODED MAGIC NUMBER. Badness.

    91    G4double pmtModuleHeight = 59.62*CLHEP::mm; //includes puck and single PMT support, not PMT base. The height of pmt module for solid works design

WCMultiPMT_container construction

There are four different versions of the container solid. The one that is used depends on two conditionals: a. mPMT_vessel_radius_curv - mPMT_vessel_radius > .1*mm b. mPMT_vessel_cyl_height < 0.01*mm

So, condition a is true and condition b is false, which puts us in the version constructed by L348-L367. The full volume is the union of a cylindrical ring and a more complicated 'cap'.

The cap is a subtraction solid in which a box volume is subtracted from a spherical shell. (This is consistent with what i see in the visualizations.) So, the first check should be whether the spherical shell of the cap has inner and outer radii consistent with those of the PMT container volume. The code:

     :
   299    G4double innerR_curv_container = mPMT_vessel_radius_curv - mPMT_outer_material_d - 54.*mm;  //-expose - dist_pmt_vessel
   300    G4double outerR_curv_container = mPMT_vessel_radius_curv - mPMT_outer_material_d;
   301
   302    //vessel_curv - vessel_height is conserved for all concentric caps !
   303    G4double innerR_container = sqrt( innerR_curv_container*innerR_curv_container -
   304                                      (mPMT_vessel_radius_curv - mPMT_vessel_cap_height)*(mPMT_vessel_radius_curv - mPMT_vessel_cap_height) );
   305    G4double outerR_container = sqrt( outerR_curv_container*outerR_curv_container -
   306                                      (mPMT_vessel_radius_curv - mPMT_vessel_cap_height)*(mPMT_vessel_radius_curv - mPMT_vessel_cap_height) );
     :
   314    G4Sphere* mPMT_top_sphere_container =
   315      new G4Sphere(    "WCmPMT_tsphere_cont",
   316                       innerR_curv_container,outerR_curv_container,
   317                       0.0*deg,360.0*deg,
   318                       0.0*deg,90.0*deg);

Ah, there is a suspect. The inner radius calculation uses a hard-coded magic number: 54.mm. In the calculation for the PMT inner radius, this term is represented by the magic number stored in the local variable pmtModuleHeight, which is 59.62mm. So, the cap shell is thinner than the pmt shell, and the PMT volumes should stick out about 5.62mm. Note that one of the reported overlaps is 5.62mm.

spradlin commented 1 year ago

Patching this problem

WCSimDetectorConstruction::ConstructPMT() and WCSimDetectorConstruction::ConstructMultiPMT() must use the same value for the shell thickness. To ensure this, they should get the value from the same source.

For now, i will make it a class constant of WCSimDetectorConstruction. This is a single-purpose solution. A more general solution would make this parameter part of the data associated with a PMT type. This probably deserves a spinoff issue.

After implementing the class constant for pmtModuleHeight, there are no more detected overlaps in the HyperK_HybridmPMT geometry, and the PMT plugs no longer extend beyond the container volume when visualized:

WCPMT_container_zoom_fix_0001

spradlin commented 1 year ago

HyperK_HybridmPMT_WithOD overlap check

Starting from my local branch with the fixes identified above for the SuperK and HyperK_HybridmPMT overlaps, run 1 event with overlap checking for the HyperK_HybridmPMT_WithOD geometry.

There are no overlaps reported by GEANT4. Great!.

spradlin commented 1 year ago

nuPRISM overlap check

Starting from my local branch with the fixes identified above for the SuperK and HyperK_HybridmPMT overlaps, run 1 event with overlap checking for the default nuPRISM geometry.

There are two overlaps detected. Here are the exceptions in the stdout:

Debug B.Q : place border pmt, spacing horiz = 285.721, vertical = -107.141
Add PMT on barrel cell 0, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (3908.4,42.6637,-150.294), overlapping by at least: 6.22355 mm 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------
Add PMT on barrel cell 0, 0
Checking overlaps for volume WCMultiPMT ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCMultiPMT
          with its mother volume WCBarrelBorderCell
          at mother local point (3921.83,-75.7705,158.186), overlapping by at least: 1.41153 cm
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

Hmm, more problems with the WCBarrelBorderCell.

spradlin commented 1 year ago

Visualization of nuPRISM WCBarrelBorderCell

nuPRISM_WCBarrelBorderCell_0000

nuPRISM_WCBarrelBorderCell_0000

I guessed this might be it. Because of the bevel, i do not think that the PMT volume will fit in the cell.

Shouldn't this have been picked up in the modulo arithmetic that set the number of PMTs in the border cell?

spradlin commented 1 year ago

Oh, and looking at this section again, i somehow completely missed the fact that the loop variables are G4doubles despite being incremented with postfix operators.

And there are some class data members that are treated as integers in context but are declared as floats.

What is this?

This can't be good.

spradlin commented 1 year ago

ToDo:

spradlin commented 1 year ago

It looks like issues with the barrel border ring geometry have been discovered and discussed before! See https://github.com/WCSim/WCSim/pull/137. Brought to my attention by @tdealtry.

spradlin commented 1 year ago

End cap construction

The 2016 merge request https://github.com/WCSim/WCSim/pull/137 very neatly describes the problems that i have been seeing with the BarrelBorderCell overlaps, and it includes a proposed solution that looks better and more durable than the patches that i have tried. The problem and its solution are summarized in the diagrams in the PDF file https://github.com/WCSim/WCSim/files/60297/OverlapBug.pdf that is attached to https://github.com/WCSim/WCSim/pull/137. Below is a JPEG copy of that diagram that should disply inline here:

OverlapBug

spradlin commented 1 year ago

End cap construction (continued)

Relation to what i have done

My solution, as described in a comment above, was to change the inter-PMT distance for the border ring so that the PMTs would fit within the border ring. This assumes that the computation of the length of the barrel border ring is what was intended.

However, after reading https://github.com/WCSim/WCSim/pull/137 and rereading relevant sections of the pack-in doc "Simulating Water Cherenkov Detectors Using WCSim" https://github.com/tdealtry/WCSim/blob/main_into_hybrid/doc/DetectorDocumentation.pdf i am convinced that this was the wrong choice. I think that the better expression of the design intention is the 'NEW' configuration in the diagram. That is, rather than changing the PMT spacing to fit in the ring, the dimensions of the ring should be set so that there is enough space for the PMTs to have the same inter-PMT spacing as in the rest of the barrel. On its interior face, the BarrelBorderCell should have at least the same 'height' as that of a standard BarrelCell (there is an additional WCBarrelPMTOffset that will be discussed below).

Quotations from "Simulating Water Cherenkov Detectors Using WCSim"

Just copying here some relevant bits of the pack-in doc, with some contextual commentary.

From Sec. 2.1 "Hierarchy of Volumes", the descriptions of the WCBarrelBorderRing and WCBarrelBorderCell volumes:

WCBarrelBorderRing This volume connects the annulus to the cap. It is essentially the uppermost (or lowermost) ring of PMTs, but has a modified bounding-box geometry and is contained in the CapAssembly because it must mate at the corner with the caps. (N.B. See 3.3, for information on corner geometry.) The border ring is divided into cells: WCBarrelBorderCell that containt he same number of PMTs (WCPMT) and the same blacksheet (WCBarrelCellBS) as the normal annulus cells.

This seems to make plain the intention that, as far as instrumentation is concerned, the WCBarrelBorderRing should be no different than a standard WCBarrelRing. It implies that whatever differences there are in the geometry are there to properly mate at the corner with the end cap.

From Sec. 3.1 "Parameters", the stated definitions of various parameters:

WCIDDiameter, WCIDHeight These two variables are used to setup the size of the detector. The height is the distance between the inner surfaces of the top and bottom blacksheets and the diameter is two times the shortest distance between the inner surface of the wall blacksheet and the center of the detector (see figure 4). This shortest radius occurs at the center of normal (not extra) cells, and is perpendicular to the blacksheet.

Specifically, the intended definition of WCIDHeight may be important---the distance between the interior surfaces of the blacksheets.

WCBarrelPMTOffset The cap volumes contain vertical space in stripes along the edges of the detector to make room for the cap PMTs (N.B. See 3.3, for information on corner geometry). This variable defines the width of these stripes. Specifically, the offset is the vertical distance from the inner surface of the cap blacksheet to the upper edge of the top cell. This edge is half the vertical PMT spacing vertically above the center of each PMT in the uppermost (or lowermost) ring.

Note that the two WCBarrelPMTOffset 'stripes' are supposed to be part of the WCIDHeight interior height. This is reinforced by the displaymath calculation in Sec. 3.2 "Example" of the document. Note also that this dimension appears as part of the height of the endcap plug in the diagram in the previous comment.

Finally, the entire text of Sec. 3.3 "Warnings" treats with exactly the kind of border problems that see,

If the PMTs have a large WCPMTExposeHeight and there is not enough space between the PMTs and the borders of the cells, the PMTs could intersect the edge of the border cells, because the border of these cells are slanted (see figure 6).

This can also happen at the caps if WCCapEdgeLimit is close to the inner radius of the detector and WCBarrelPMTOffset is small.

During the setup an incomplete checking for obvious overlaps occurs. You should see a lot of the flowing lines:
Checking overlaps for volume WCBarrelPMT ... OK!
If you see warnings instead, it is likely that there are too large or too many PMTs. An absence of warnings does not mean that there are no overlaps.

There is no check if the placement of the PMTs on caps is correct. It would take to much time to do this check every time, because there are too many PMTs in side a single volume.

Some rules of thumb to avoid overlaps are: WCBarrelPMTOffset > WCPMTExposeHeight or WCIDDiameter / 2 - WCCapEdgeLimit > WCPMTExposeHeight ensures cap PMTs are fullywithinthecapvolume,andvertical spacing of the PMTs / 2 > WCPMTExposeHeight + WCPMTRadius ensures the top and bottom ring PMTs are within WCBarrelBorderRing. These rules are general, and PMTs may be placed closer with careful attention to geometry. Collisions are not an issue in Super-K, and newer detectors with smaller, more-efficient (and thus sparser) PMTs should have little issue provided reasonable specifications of WCCapEdgeLimit and WCBarrelPMTOffset.

spradlin commented 1 year ago

Computation of ring heights

The computation of the height of standard barrel cell is in WCSimConstructCylinder.cc L89:

  barrelCellHeight  = (WCIDHeight-2.*WCBarrelPMTOffset)/WCBarrelNRings;

The calculation clearly indicates that the WCIDHeight includes WCBarrelNRings instrumented rings and the two buffer 'stripes' at the ends of the barrel that have a height of WCBarrelPMTOffset each.

The computation of the derived mainAnnulusHeight in WCSimConstructCylinder.cc L91, which becomes the height of the central section of the barrel, is

  mainAnnulusHeight = WCIDHeight -2.*WCBarrelPMTOffset -2.*barrelCellHeight;

which also indicates that the interior side of the border ring should be (barrelCellHeight + WCBarrelPMTOffset).

OK, so far. Now let me go back into the border ring assembly in WCSimDetectorConstruction::ConstructCaps(). The dimensions are defined in parallel arrays that are described in a previous comment and reproduced here:

    39  #define MIRROR_WCSIM_DEVELOP_POLY
     :
  1447    G4double borderAnnulusZ[3] = {
  1448  #ifdef MIRROR_WCSIM_DEVELOP_POLY
  1449          -barrelCellHeight/2.*zflip,
  1450          (-barrelCellHeight/2.+(WCIDRadius-innerAnnulusRadius))*zflip,
  1451  #else
  1452          (-barrelCellHeight/2.-(WCIDRadius-innerAnnulusRadius))*zflip,
  1453          -barrelCellHeight/2.*zflip,
  1454  #endif
  1455          barrelCellHeight/2.*zflip};
  1456    G4double borderAnnulusRmin[3] = { WCIDRadius, innerAnnulusRadius, innerAnnulusRadius};
  1457    G4double borderAnnulusRmax[3] = {outerAnnulusRadius, outerAnnulusRadius,outerAnnulusRadius};

Here, the height of the interior face is either barrelCellHeight or (barrelCellHeight -WCIDRadius - innerAnnulusRadius ) rather than (barrelCellHeight + WCBarrelPMTOffset). The WCBarrelPMTOffset seems to be incorporated into the cap part. Ooh, this may take some disentangling.

Or, perhaps not. If the WCBarrelPMTOffset stripe is incorporated into the cap volume, and if the inner height of the barrel border ring is already equal to barrelCellHeight, then it may be that WCBarrelPMTOffset can be removed from all of the calculations related to the barrel border ring and produce the correct results. I may need to undo some of the changes that i initially made to address the problem.

Before addressing that, what is going on with this conditionally compiled code? The #else clause gives the interior face height barrelCellHeight. The #ifdef clause carves the bevel out of the cell height so that the total height is barrelCellHeight. So, turning off MIRROR_WCSIM_DEVELOP_POLY seems to indicated. But, what is this condition supposed to accomplish? This is followed up in a subsequent comment.

spradlin commented 1 year ago

Revisiting the WCBorderPMTOffset

The solution that i proposed in a previous comment was to update the value of WCSimDetectorConstruction::WCBorderPMTOffset to make it large enough to exclude the bevel of the border cell. Interestingly, i have been conflating WCSimDetectorConstruction::WCBorderPMTOffset and WCSimDetectorConstruction::WCBarrelPMTOffset in my thinking to this point. However, they are distinct variables, and they have independent uses.

While WCSimDetectorConstruction::WCBarrelPMTOffset has a defined meaning in the documentation, WCBorderPMTOffset does not. So, what is this WCBorderPMTOffset? Where did it come from? Why is it not just the height of the bevel?

Interestingly, it is only defined in two detector geometry configurations:

   363  void WCSimDetectorConstruction::SetHyperKWithODGeometry()
   364  {
   365    WCDetectorName = "HyperKWithOD";
    :
   374    WCBarrelPMTOffset     = WCPMTRadius; //offset from vertical
   375    WCBorderPMTOffset     = 0.25*m; // PMT offset in border cells. Also makes distance between PMTs shorter.
    :
   443  }
    :
   752  void WCSimDetectorConstruction::SetHyperK_HybridmPMT_WithOD_Geometry()
   753  {
   754    //Box&Line
   755    WCDetectorName = "HyperK_HybridmPMT_WithOD";
    :
   798    WCBarrelPMTOffset     = std::max(WCPMTRadius,mPMT_vessel_tot_height); //offset from vertical
   799    WCBorderPMTOffset     = 0.25*m; // PMT offset in border cells. Also makes distance between PMTs shorter.
    :
   875  } 

The comments associated to the WCBorderPMTOffset assignments are interesting, and contrary to the understanding that i have come to regarding the instrumentation of the border rings. Where did this come from?

Using GiHub's blame mechanism, i tracked the introduction of to this commit: https://github.com/tdealtry/WCSim/commit/be3bc48e6e8ff752abe1494deffc3966b13897d8. The comment attached to the commit expands on the comments in the code:

Added a variable to offset PMTs in the border ring. This stops them from overlapping with the slanted angle. The offset also reduces the overall height of the cell so PMTs are slightly closer together so that they can all fit.

Ah, so it was an incomplete implementation of a mechanism to solve the problem. I am convinced, now, that this introduces additional complications that should be unnecessary.

spradlin commented 1 year ago

MIRROR_WCSIM_DEVELOP_POLY conditionals

Use of the compiler macro MIRROR_WCSIM_DEVELOP_POLY is confined to src/WCSimConstructCylinder.cc:

$ grep --with-filename --line-number MIRROR_WCSIM_DEVELOP_POLY include/*.hh src/*.cc
src/WCSimConstructCylinder.cc:39:#define MIRROR_WCSIM_DEVELOP_POLY
src/WCSimConstructCylinder.cc:1448:#ifdef MIRROR_WCSIM_DEVELOP_POLY
src/WCSimConstructCylinder.cc:1648:#ifdef MIRROR_WCSIM_DEVELOP_POLY
src/WCSimConstructCylinder.cc:1768:#ifdef MIRROR_WCSIM_DEVELOP_POLY

So, let us take a look at the three segments of code that test the macro.

The first is the definition of the border ring dimensions included in a previous comment.

The second is part of the construction of the WCCap volume in WCSimConstructCylinder.cc L1643-L1715. The WCCap logical volume is the volume in which the endcap PMTs and the endcap blacksheet are embedded. The dimensions of the cap volume are defined in arrays that are intended to mate with the WCBarrelBorderRing:

  1647    G4double capZ[4] =
  1648  #ifdef MIRROR_WCSIM_DEVELOP_POLY
  1649          { (-WCBlackSheetThickness-1.*mm)*zflip,
  1650            WCBarrelPMTOffset*zflip,
  1651            WCBarrelPMTOffset*zflip,
  1652            (WCBarrelPMTOffset+(WCIDRadius-innerAnnulusRadius))*zflip} ;
  1653  #else   
  1654          { (-WCBlackSheetThickness-1.*mm)*zflip,
  1655            (WCBarrelPMTOffset - (WCIDRadius-innerAnnulusRadius))*zflip,
  1656            (WCBarrelPMTOffset - (WCIDRadius-innerAnnulusRadius))*zflip,
  1657            WCBarrelPMTOffset*zflip} ;
  1658  #endif
  1659    
  1660    G4double capRmin[4] = {  0. , 0., 0., 0.} ;
  1661    G4double capRmax[4] = {outerAnnulusRadius, outerAnnulusRadius,  WCIDRadius, innerAnnulusRadius};

The third instance is part of the construction of the WCCapBlackSheet in WCSimConstructCylinder.cc L1767-L1822. This makes sense, as the blacksheet dimensions need to match. Like the WCBarrelBorderRing and WCCap, the WCCapBlackSheet is defined as a G4Polyhedra with dimensions controlled by arrays:

  1767    G4double capBlackSheetZ[4] = {-WCBlackSheetThickness*zflip, 0., 0.,
  1768  #ifdef MIRROR_WCSIM_DEVELOP_POLY
  1769                                                                  WCBarrelPMTOffset*zflip}; 
  1770  #else     
  1771                                                                  (WCBarrelPMTOffset - (WCIDRadius-innerAnnulusRadius)) *zflip};
  1772  #endif    
  1773    G4double capBlackSheetRmin[4] = {0., 0., WCIDRadius, WCIDRadius};
  1774    G4double capBlackSheetRmax[4] = {WCIDRadius+WCBlackSheetThickness,
  1775                                     WCIDRadius+WCBlackSheetThickness,
  1776                                                                     WCIDRadius+WCBlackSheetThickness,
  1777                                                                     WCIDRadius+WCBlackSheetThickness};

Note that no geometry elements are ever constructed at part of the conditionally compiled code. The dimensions are changed for no apparent reason. OK, let's excise it. Will probably need to recheck the alignment of everything.

spradlin commented 1 year ago

Removing WCBorderPMTOffset

Right, so let us get rid of it. Certainly, eliminate it as a data member of WCSimDetectorConstruction. It does have a place as a local variable within WCSimDetectorConstruction::ConstructCaps(), but its value should always be the height of the bevel.

Further, this change to the PMT vertical spacing should not be changed---it should be identical to that of a normal border cell. The offset should not figure into the computation of the spacing. The offset should just be exactly the height of the bevel, which means that the PMT placement in the border cell should fit as well as they do in a standard barrel cell.

The offset should only appear in the computation of the placement vector.

spradlin commented 1 year ago

ToDo:

More things to check:

spradlin commented 1 year ago

Blacksheet in WCBarrelBorderCell

After the removal of the compiler conditional and the removal of WCBorderPMTOffset, the blacksheet volumes no longer cover the full height of the cell, as shown in the visualization below for the WCBarrelBorderCell of the HyperK_HybridmPMT_WithOD geometry:

HyperK_HybridmPMT_WithOD_WCBarrelBorderCell_0000

The green box outlines near the perimeter of the cell indicate the extent of the blacksheets---both the inner detector blacksheet and outer detector Tyvek are shown. As can be seen, they stop where the blacksheets for a normal BarrelCell would stop. They do not extend to the edge of the bevel.

This is simply because copies of the normal BarrelCell blacksheet and normal BarrelCell OD Tyvek are placed in these cells.

Custom blacksheets must be constructed for the border ring cells.

Right, so i can copy-mod the code for normal BarrelCells from WCSimDetectorConstruction::ConstructCylinder() that constructs these volumes. A better method might be to make a cell constructing function or class that can build these things self-consistently without code duplication. I think that would involve a substantial redesign of the geometry construction code---although it is probably worth doing at some point. For now, let us go with the hacky copy-mod solution. Need to do this four times: (ID, OD) x (standard, extra tower). Ugh.

Ok, so after fixing the blacksheets for the border cells, i get the following updated visualization of a WCBarrelBorderCell of HyperK_HybridmPMT_WithOD:

HyperK_HybridmPMT_WithOD_WCBarrelBorderCell_after

spradlin commented 1 year ago

Visualizations of WCBarrelBorderCell in other geometries

After the changes to the blacksheet dimensions described in the previous comment, produce visualizations of other geometries.

SuperK

SuperK_WCBarrelBorderCell_after

HyperK_HybridmPMT

HyperK_HybridmPMT_WCBarrelBorderCell_after

nuPRISM

nuPRISM_WCBarrelBorderCell_after

The overlaps that prompted this line of modifications are no longer present.

Good.

spradlin commented 1 year ago

Visualizations of WCExtraTowerBorderCell

Check the extra cell in the border ring for geometries that have one.

SuperK

SuperK_WCExtraTowerBorderCell_0000

HyperK_HybridmPMT

HyperK_HybridmPMT_WCExtraTowerBorderCell_0000

HyperK_HybridmPMT_WithOD

HyperK_HybridmPMT_WithOD_WCExtraTowerBorderCell_0000

That looks good.

spradlin commented 1 year ago

Checking OD PMT placement in WCBarrelBorderCell

Like the ID PMT placement, the OD PMT placement should be the same as for a regular BarrelCell. The placement for regular BarrelCells is in WCSimConstructCylinder.cc L1201-L1228:

  1167      G4double barrelODCellWidth   = 2.*WCODRadius*tan(dPhi/2.);
  1168      G4double barrelODCellHeight  = barrelCellHeight * (barrelODCellWidth/barrelCellWidth);
   :
  1201      G4double horizontalODSpacing = barrelODCellWidth/WCPMTODperCellHorizontal;  
  1202      G4double verticalODSpacing   = barrelODCellHeight/WCPMTODperCellVertical;
  1203  
  1204      if(WCODPMTShift > barrelODCellWidth/2. - WCPMTODRadius) WCODPMTShift = 0.*cm;
  1205  
  1206      for(G4long i = 0; i < WCPMTODperCellHorizontal; i++){
  1207        for(G4long j = 0; j < WCPMTODperCellVertical; j++){
  1208                  G4cout << "Adding OD PMT in cell " << i << ", " << j << G4endl; 
  1209          G4ThreeVector Container =  G4ThreeVector(WCODRadius,
  1210                                                   -barrelODCellWidth/2.+(i+0.5)*horizontalODSpacing+((G4int)(std::pow(-1,j))*(G4int)(WCODPMTShift)/2),
  1211                                                   -(barrelCellHeight * (barrelODCellWidth/barrelCellWidth))/2.+(j+0.5)*verticalODSpacing);
  1212
  1213                  //              G4cout << " qqqqqqqqqqqqqqqqqqqqqqqq barrel i " << i << " of " << WCPMTODperCellHorizontal << " j " << j << " of " << WCPMTODperCellVertical << " Container (" << Container.x() << ", " << Container.y()
  1214                  //                                << ", " << Container.z() << ") " << G4endl;
  1215
  1216          //G4VPhysicalVolume* physiWCBarrelWLSPlate =
  1217                  new G4PVPlacement(WCPMTODRotation,           // its rotation
  1218                                                    Container,
  1219                                                    logicWCODWLSAndPMT,         // its logical volume
  1220                                                    "WCBarrelCellODContainer",  // its name
  1221                                                    logicWCBarrelCell,         // its mother volume
  1222                                                    false,                     // no boolean operations
  1223                                                    i*WCPMTODperCellVertical+j,
  1224                                                    checkOverlapsPMT);
  1225
  1226
  1227        }
  1228      }

The OD PMT placement for BarrelBorderCells is in WCSimConstructCylinder.cc L2119-L2144:

  2119      const G4double barrelODCellWidth   = 2.*WCODRadius*tan(dPhi/2.);
  2120      const G4double barrelODCellHeight  = barrelCellHeight * (barrelODCellWidth/barrelCellWidth);
  2121      G4double horizontalODSpacing = barrelODCellWidth/WCPMTODperCellHorizontal;
  2122      const G4double verticalODSpacing   = barrelODCellHeight / WCPMTODperCellVertical;             
  2123                          
  2124      for(G4long i = 0; i < WCPMTODperCellHorizontal; i++){
  2125        for(G4long j = 0; j < WCPMTODperCellVertical; j++){ 
  2126          G4ThreeVector Container =  G4ThreeVector(WCODRadius,
  2127                                                   -barrelODCellWidth/2.+(i+0.5)*horizontalODSpacing+((G4int)(std::pow(-1,j))*(G4int)(WCODPMTShift)/2),
  2128                                                   -(barrelCellHeight * (barrelODCellWidth/barrelCellWidth))/2.+(j+0.5)*verticalODSpacing); 
  2129  
  2130                  G4cout << "Adding OD PMT in barrel in cell" << i << ", " << j << G4endl;
  2131          //G4VPhysicalVolume* physiWCBarrelWLSPlate =      
  2132                  new G4PVPlacement(WCPMTODRotation,           // its rotation
  2133                                                    Container,
  2134                                                    logicWCODWLSAndPMT,         // its logical volume
  2135                                                    "WCBorderCellODContainer",  // its name
  2136                                                    logicWCBarrelBorderCell,         // its mother volume                                   
  2137                                                    false,                     // no boolean operations                                     
  2138                                                    i*WCPMTODperCellVertical+j,
  2139                                                    checkOverlapsPMT);
  2140                                                            
  2141           
  2142                                                            
  2143        }                                                   
  2144      }

That looks consistent. There is no zflip in the placement vector in the logicWCBarrelBorderCell. Should there be? I believe so---to check.

Aaaargh! The reuse of names makes it hard to check components of the top and bottom caps independently. I may need to disambiguate the names for the top and bottom cap components.

It seems like GEANT4 would have a reflection transform that would make the duplicate construction unecessary....

spradlin commented 1 year ago

Checking PMT placement in the extra tower

Again, the PMT placement in the extra tower cells should be uniform throughout the barrel. Placement of the ID PMTs in the ExtraTowerCell is in WCSimConstructCylinder.cc L1006-L1035:

  1006            G4double towerWidth = WCIDRadius*tan(2*pi-totalAngle);
  1007
  1008            G4double horizontalSpacingExtra   = towerWidth/(WCBarrelNumPMTHorizontal-WCBarrelRingNPhi*WCPMTperCellHorizontal);
  1009            // verticalSpacing is identical to that for normal cell.
  1010
  1011            for(G4long i = 0; i < (WCBarrelNumPMTHorizontal-WCBarrelRingNPhi*WCPMTperCellHorizontal); i++){
  1012                  for(G4long j = 0; j < WCPMTperCellVertical; j++){
  1013                    G4ThreeVector PMTPosition =  G4ThreeVector(WCIDRadius/cos(dPhi/2.)*cos((2.*pi-totalAngle)/2.),
  1014                                                                                                           towerWidth/2.-(i+0.5)*horizontalSpacingExtra,
  1015                                                                                                           -barrelCellHeight/2.+(j+0.5)*verticalSpacing);
  1016                    PMTPosition.rotateZ(-(2*pi-totalAngle)/2.); // align with the symmetry                          
  1017                    //axes of the cell 
  1018      
  1019                    //G4cout << "Adding PMT in extra tower " << i << ", " << j << " position: " << PMTPosition << G4endl;           
  1020  #ifdef ACTIVATE_IDPMTS                                                  
  1021                    //G4VPhysicalVolume* physiWCBarrelPMT =
  1022                    new G4PVPlacement(WCExtraPMTRotation,             // its rotation                                               
  1023                                                          PMTPosition,    
  1024                                                          (i==j && hybrid && WCPMTPercentCoverage2!=0)?logicWCPMT2:logicWCPMT,                // its logical volume
  1025                                                          pmtname,             // its name
  1026                                                          logicWCExtraTowerCell,         // its mother volume
  1027                                                          false,                     // no boolean operations
  1028                                                          i*WCPMTperCellVertical+j,
  1029                                                          checkOverlapsPMT);
  1030  #endif
  1031                    // logicWCPMT->GetDaughter(0),physiCapPMT is the glass face. If you add more
  1032                    // daughter volumes to the PMTs (e.g. a acryl cover) you have to check, if
  1033                    // this is still the case.
  1034                  }
  1035            }

In the border ring, ID PMT placement in ... is WCSimConstructCylinder.cc L2045-L2076:

  2045            G4double towerWidth = WCIDRadius*tan(2*pi-totalAngle);
  2046
  2047            G4double horizontalSpacingExtra   = towerWidth/(WCBarrelNumPMTHorizontal-WCBarrelRingNPhi*WCPMTperCellHorizontal);
  2048            // verticalSpacing is identical to that for normal cell.
  2049                    
  2050            for(G4long i = 0; i < (WCBarrelNumPMTHorizontal-WCBarrelRingNPhi*WCPMTperCellHorizontal); i++){
  2051                  for(G4long j = 0; j < WCPMTperCellVertical; j++){
  2052                    G4ThreeVector PMTPosition =  G4ThreeVector(WCIDRadius/cos(dPhi/2.)*cos((2.*pi-totalAngle)/2.),
  2053                                                                                                           towerWidth/2.-(i+0.5)*horizontalSpacingExtra,  
  2054                                                                                                           (-(barrelCellHeight-WCBorderPMTOffset)/2.+(j+0.5)*verticalSpacing)*zflip);
  2055                    PMTPosition.rotateZ(-(2*pi-totalAngle)/2.); // align with the symmetry 
  2056                    //axes of the cell                    
  2057  #ifdef ACTIVATE_IDPMTS
  2058  #ifdef WCSIMCONSTRUCTCYLINDER_VERBOSE                   
  2059                    G4cout << "Add PMTs in extra tower, cell " << i << ", " << j << G4endl;
  2060  #endif
  2061                    //G4VPhysicalVolume* physiWCBarrelBorderPMT =
  2062                    new G4PVPlacement(WCExtraPMTRotation,                          // its rotation
  2063                                                          PMTPosition,
  2064                                                          (i==j && hybrid && WCPMTPercentCoverage2!=0)?logicWCPMT2:logicWCPMT,                 // its logical volume
  2065                                                          "WCPMT",             // its name
  2066                                                          logicWCExtraBorderCell,         // its mother volume
  2067                                                          false,                     // no boolean operations
  2068                                                          i*WCPMTperCellVertical+j,
  2069                                                          checkOverlapsPMT);
  2070  #endif
  2071
  2072                    // logicWCPMT->GetDaughter(0),physiCapPMT is the glass face. If you add more
  2073                    // daugter volumes to the PMTs (e.g. a acryl cover) you have to check, if
  2074                    // this is still the case.
  2075                  }
  2076            }

The ID placement looks consistent.

spradlin commented 1 year ago

Endcap component names

The implementation of zflip in WCSimDetectorConstruction::ConstructCaps() is complicated. In order to visually debug it, i need to be able to visualize the components of the top and bottom cap independently. This means that corresponding components in the top and bottom caps need different names. Let me collect the names of geometry elements constructed in WCSimDetectorConstruction::ConstructCaps() and see what modifications might make sense to indicate chirality.

How about updating substrings

spradlin commented 1 year ago

Chiral visualizations: SuperK

After updating the names of the geometry elements in the endcaps to distinguish parts constructed for the Top and Bottom endcaps, i can make visualizations of the cells of both border rings.

WCBarrelTopBorderCell

SuperK_WCBarrelTopBorderCell_0000

WCBarrelBotBorderCell

SuperK_WCBarrelBotBorderCell_0000

WCExtraTowerTopBorderCell

SuperK_WCExtraTowerTopBorderCell_0000

WCExtraTowerBotBorderCell

SuperK_WCExtraTowerBotBorderCell_0000

spradlin commented 1 year ago

Chiral visualizations: HyperK_HybridmPMT_WithOD

WCBarrelTopBorderCell

HyperK_HybridmPMT_WithOD_WCBarrelTopBorderCell_0000

WCBarrelBotBorderCell

HyperK_HybridmPMT_WithOD_WCBarrelBotBorderCell_0000

WCExtraTowerTopBorderCell

HyperK_HybridmPMT_WithOD_WCExtraTowerTopBorderCell_0000

WCExtraTowerBotBorderCell

HyperK_HybridmPMT_WithOD_WCExtraTowerBotBorderCell_0000

spradlin commented 1 year ago

Updated overlap checks

This is all looking good. Let us check for overlaps again. With overlap checking turned on, run a 1 event job for each of several geometries. Check the stdout/stderr for detected overlaps.

The overlap reports for HyperK_HybridmPMT are

Checking overlaps for volume WCTopCap ... OK! 
Checking overlaps for volume WCTopCapBlackSheet ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCTopCapBlackSheet
          with its mother volume WCTopCapPolygon 
          at mother local point (18098.5,26951.8,0.22657), overlapping by at least: 546.898 um 
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------
Checking overlaps for volume WCBotCap ... OK! 
Checking overlaps for volume WCBotCapBlackSheet ...
-------- WWWW ------- G4Exception-START -------- WWWW -------
*** G4Exception : GeomVol1002
      issued by : G4PVPlacement::CheckOverlaps()
Overlap with mother volume !
          Overlap is detected for volume WCBotCapBlackSheet
          with its mother volume WCBotCapPolygon 
          at mother local point (-32005,5045.06,-0.936551), overlapping by at least: 44.8654 um
NOTE: Reached maximum fixed number -1- of overlaps reports for this volume !
*** This is just a warning message. ***
-------- WWWW -------- G4Exception-END --------- WWWW -------

The overlaps are small, but it looks like i may have missed something when the height of the border ring changed. Let us look at the cap black sheet construction and placement. A first set of suspects are the conditionally compiled dimensions in a previous comment, where i have removed the conditionals.

spradlin commented 1 year ago

More end cap construction

Now that i am convinced that the bulk of the barrel border ring is constructed as intended, let us revisit the remainder of the construction and make sure that things are fitting appropriately. I will rely heavily on the diagram attached to https://github.com/WCSim/WCSim/pull/137 and reproduced in a previous comment. Specifically, i will focus on the height alignment.

Containing volume

The entire endcap assembly, including the border ring, is included in a cylindrical volume called the [Top|Bot]CapAssembly (CapAssembly before my name changes):

  1425    capAssemblyHeight = (WCIDHeight-mainAnnulusHeight)/2+1*mm+WCBlackSheetThickness;
  1426  
  1427    G4Tubs* solidCapAssembly = new G4Tubs("CapAssembly",
  1428                                                                                  0.0*m,
  1429                                                                                  outerAnnulusRadius/cos(dPhi/2.),
  1430                                                                                  capAssemblyHeight/2,
  1431                                                                                  0.*deg,
  1432                                                                                  360.*deg);

WCIDHeight and WCBlackSheetThickness are basic configuration parameters; the mainAnnulusHeight is a computed value:

    91    mainAnnulusHeight = WCIDHeight -2.*WCBarrelPMTOffset -2.*barrelCellHeight;

so
(WCIDHeight-mainAnnulusHeight)/2 = WCBarrelPMTOffset + barrelCellHeight,
and the total height of the endcap containing volume should be
capAssemblyHeight = barrelCellHeight + WCBarrelPMTOffset + WCBlackSheetThickness + 1mm
which is consistent with the diagram with a 1mm buffer.

Border ring

The first subvolume constructed and placed in the CapAssembly is the barrel border ring, which was my focus up to now. After my fixes to the barrel border ring, its height profile is computed as

  G4double borderAnnulusZ[3] = {
        (-barrelCellHeight/2.-(WCIDRadius-innerAnnulusRadius))*zflip,
        -barrelCellHeight/2.*zflip,
        barrelCellHeight/2.*zflip};

so, having a height of
barrelCellHeight
on its interior side, but increasing to
barrelCellHeight + (WCIDRadius-innerAnnulusRadius)
with the addition of the bevel.

OK, so there will be problems if the ring ends up taller than the cap assembly, that is, if
(WCIDRadius-innerAnnulusRadius) > WCBarrelPMTOffset + WCBlackSheetThickness + 1mm
or, because one WCBlackSheetThickness is needed above the ring, there will be problems if
(WCIDRadius-innerAnnulusRadius) > WCBarrelPMTOffset + 1mm
Are there conditions under which this might happen? I have looked at innerAnnulusRadius in a previous comment. I did some algebra to determine that
WCIDRadius - innerAnnulusRadius = WCPMTExposeHeight + 1mm
so we must have
WCBarrelPMTOffset + WCBlackSheetThickness >= WCPMTExposeHeight

Hmm, that is not precisely the case for the HyperK_HybridmPMT geometries. For those, we have the following conditional definition of innerAnnulusRadius

   100    if(hybrid){
   101          if(mPMT_vessel_cyl_height + mPMT_vessel_radius < 1.*mm)
   102            innerAnnulusRadius = WCIDRadius - std::max(WCPMTExposeHeight,WCPMTExposeHeight2)-1.*mm;
   103          else
   104            innerAnnulusRadius = WCIDRadius - std::max(WCPMTExposeHeight,(mPMT_vessel_cyl_height + mPMT_vessel_radius)) -1.*mm;
   105    }

Pulling the numbers for the HyperK_HybridmPMT geometry,

   531  void WCSimDetectorConstruction::SetHyperK_HybridmPMTGeometry()
   532  {
    :
   536    WCSimPMTObject * PMT = CreatePMTObject("BoxandLine20inchHQE", WCIDCollectionName);
    :
   538    WCPMTExposeHeight   = PMT->GetExposeHeight();
   539    WCPMTRadius         = PMT->GetRadius();
    :
   545    mPMT_ID_PMT = "PMT3inchR14374";    //can be changed in macro through mPMT settings.
    :
   547    WCSimPMTObject * PMT2 = CreatePMTObject(mPMT_ID_PMT, WCIDCollectionName2);
    :
   549    WCPMTExposeHeight2 = PMT2->GetExposeHeight();
    :
   554    mPMT_vessel_cyl_height = 38.*CLHEP::mm;    //option A, option B would be 277 mm
    :
   556    mPMT_vessel_radius = 254.*CLHEP::mm;
    :
   572    G4double mPMT_vessel_tot_height = mPMT_vessel_radius + mPMT_vessel_cyl_height;
    :
   577    WCBarrelPMTOffset     = std::max(WCPMTRadius,mPMT_vessel_tot_height); //offset from vertical
    :
   595  } 

In WCSimPMTObject.cc

   :
  2115  G4double BoxandLine20inchHQE::GetExposeHeight() {return .18*m;}
  2116  G4double BoxandLine20inchHQE::GetRadius() {return .254*m;}
   :
  3152  G4double PMT3inchR14374::GetExposeHeight() {return 0.02*m;}          //.0153*m;} //.0200*m;}  //from TechSheet for 3in (only photocathode would be 15.3mm h, for a radius as photocathode of 36 mm)
   :

First, calculate WCBarrelPMTOffset

  mPMT_vessel_tot_height = mPMT_vessel_radius + mPMT_vessel_cyl_height
      = 254 mm + 38 mm
      = 292 mm
  WCBarrelPMTOffset = max(WCPMTRadius, mPMT_vessel_tot_height)
      = max(254 mm, 292 mm)
      = 292 mm

Now, mPMT_vessel_cyl_height + mPMT_vessel_radius is > 1 mm, so we are in the else clause:

  innerAnnulusRadius = WCIDRadius - max(WCPMTExposeHeight, (mPMT_vessel_cyl_height + mPMT_vessel_radius))  - 1 mm
  WCIDRadius - innerAnnulusRadius = max(WCPMTExposeHeight, (mPMT_vessel_cyl_height + mPMT_vessel_radius))  + 1 mm
      = max(180 mm, 292 mm) + 1 mm
      = 293 mm

The condition that we need to satisfy is

  WCBarrelPMTOffset + 1 mm >= WCIDRadius - innerAnnulusRadius
  292 mm + 1 mm >= 293 mm

So, neglecting the floating point precision, the equality holds. This might be ok.

It makes sense to put in an assert to enforce this. Hmm, assert is not working as i expected it to. It might pay to look into setting up some dimension checking exceptions in the geometry construction. For now, let me put in an error message if the total height of the border ring (with its bevel) exceeds the height of the CapAssembly containing volume.

This new test triggers for the HyperK_HybridmPMT geometries! The reporting that i put in gives me the following:

IMPOSSIBLE GEOMETRY:  capAssemblyHeight (4967.79) too small to contain border ring (4947.79) and endcap blacksheet (20)

I would guess that floating point precision differences are at play. Equality does not provide enough buffer for floating point errors. I suspect that this can be easily fixed by inflating WCBarrelPMTOffset a little for the HyperK_HybridmPMT geometries. Multiply it by (1 + epsilon) and see if that corrects it. It seems to.

Does slightly inflating WCBarrelPMTOffset also abrogate the overlaps? No. This is not surprising. The relative size of the overlaps is much larger than epsilon.

Continuing.

The border ring should be placed within the CapAssembly so that the face opposite that with the bevel is conincident with one of the faces of the CapAssembly. The placement code:

  1474    new G4PVPlacement(0,
  1475                                          G4ThreeVector(0.,0.,(capAssemblyHeight/2.- barrelCellHeight/2.)*zflip),
  1476                                          logicWCBarrelBorderRing,
  1477                                          "WCBarrelBorderRing",
  1478                                          logicCapAssembly,
  1479                                          false, 0,
  1480                                          checkOverlaps);

Hmm, that superficially makes sense, but i will need to check the orientation with a visualization. This would probably be a good point at which to look at endcap assembly.

Turning off the PMT placement to reduce visual noise, let's look at an endcap. Here is a section of the TopCapAssembly of the nuPRISM geometry:

nuPRISM_TopCapAssembly_0000

That looks like the border ring orientation makes sense.

Extra tower cell

The next geometry element placed in the CapAssembly is WCExtraTower[Top|Bot]BorderCell, if one is necessary. It is constructed with the same z profile as the border ring. It should be placed at the same height, as well:

  1610      physiWCExtraBorderCell =
  1611        new G4PVPlacement(0,
  1612                                                  G4ThreeVector(0.,0.,(capAssemblyHeight/2.- barrelCellHeight/2.)*zflip),
  1613                                                  logicWCExtraBorderCell,
  1614                                                  "WCExtraTowerBorderCell",
  1615                                                  logicCapAssembly,
  1616                                                  false, 0,
  1617                                                  checkOverlaps);

which it is.

The last geometry element placed in the CapAssembly is the Cap, the plug assembly. This contains the BlackSheet. Look at its construction in the next comment.

Z so far

Ok, so my understanding of the z extents so far. In the coordinate system of the CapAssembly, the containing volume extends from
z_{Cap,min} = -capAssemblyHeight/2 to z_{Cap,max} = capAssemblyHeight/2.
For now, i will just treat the zflip = +1 case. In the coordinate system of the WCBarrelBorderRing, it extends from z'_{Ring,min} = -barrelCellHeight/2.-(WCIDRadius-innerAnnulusRadius) to
z'_{Ring,max} = barrelCellHeight/2.
The transformation to the coordinate system of the CapAssembly should be a simple translation in z that is defined by the placement vector. The placement vector sets
z' = 0 at z = capAssemblyHeight/2.- barrelCellHeight/2.
so that
z = z' + capAssemblyHeight/2.- barrelCellHeight/2.
= z' + ( WCBarrelPMTOffset + WCBlackSheetThickness + 1mm)/2. Applying this transformation gives the extent of the WCBarrelBorderRing in the coordinate system of the CapAssembly as
z_{Ring,min} = capAssemblyHeight/2. - barrelCellHeight - (WCIDRadius-innerAnnulusRadius) to
z_{Ring,max} = capAssemblyHeight/2.
The end of the bevel of the WCBarrelBorderRing is at a z position of
z_{Ring,bev} = capAssemblyHeight/2. - barrelCellHeight
exactly one barrelCellHeight below the face of the CapAssembly. Good, i think.

spradlin commented 1 year ago

More end cap construction (cont'd)

The Cap (plug volume)

The cap volume is a plug that contains the endcap PMTs and the blacksheet behind them. How it is constructed depends on whether there is an ExtraTower in the barrel:

The z profile of WC[Top|Bot]Cap is the same for both cases. After removing conditional code, the definition of the z profile is

  G4double capZ[4] =
        { (-WCBlackSheetThickness-1.*mm)*zflip,
          (WCBarrelPMTOffset - (WCIDRadius-innerAnnulusRadius))*zflip,
          (WCBarrelPMTOffset - (WCIDRadius-innerAnnulusRadius))*zflip,
          WCBarrelPMTOffset*zflip} ;
  G4double capRmin[4] = {  0. , 0., 0., 0.} ;
  G4double capRmax[4] = {outerAnnulusRadius, outerAnnulusRadius,  WCIDRadius, innerAnnulusRadius};

Its placement in the CapAssembly is also independent of its geometric construction:

  1723    G4VPhysicalVolume* physiWCCap =               
  1724      new G4PVPlacement(0,                           // no rotation
  1725                                            G4ThreeVector(0.,0.,(-capAssemblyHeight/2.+1*mm+WCBlackSheetThickness)*zflip),     // its position
  1726                                            logicWCCap,          // its logical volume    
  1727                                            "WCCap",             // its name
  1728                                            logicCapAssembly,                  // its mother volume
  1729                                            false,                       // no boolean operations
  1730                                            0,                          // Copy #                                                 
  1731                                            checkOverlaps);

z in CapAssembly

Let's run the arithmetic. In the coordinate system of the WC[Top|Bot]Cap plug volume, the extent in z'' for zflip = +1 is
z''_{plug,min} = -WCBlackSheetThickness - 1.*mm
z''_{plug,max} = WCBarrelPMTOffset
with an intermediate position for the bevel edge of
z''_{plug,bev} = WCBarrelPMTOffset - (WCIDRadius - innerAnnulusRadius)

The transformation to the CapAssembly coordinate system implied by its placement is
z = z'' - capAssemblyHeight/2. + 1*mm + WCBlackSheetThickness
so that the plug in CapAssembly coordinates extends
z_{plug,min} = -capAssemblyHeight/2.
z_{plug,max} = -capAssemblyHeight/2. + WCBarrelPMTOffset + WCBlackSheetThickness + 1*mm
The minimum extent corresponds to a face of the CapAssembly. Good, i think.

In order for the plug to align properly with the border ring, we should expect that
z_{plug,max} == z_{Ring,bev} and z_{plug,bev} == z_{Ring,min}.
So, check this.

  z_{plug,max} =? z_{Ring,bev}
  -capAssemblyHeight/2. + WCBarrelPMTOffset + WCBlackSheetThickness + 1*mm =? capAssemblyHeight/2. - barrelCellHeight
  capAssemblyHeight =? barrelCellHeight + WCBarrelPMTOffset + WCBlackSheetThickness + 1*mm

Yes! Good. That looks correct.

End black sheet

The CapBlackSheet should be, approximately, a hollowed-out copy of the plug volume that stops at the bevel. It construction is dependent on whether an ExtraTower is constructed, like the plug. However, its z profile is independent of this consideration. After removing conditional code, the definition of the z profile is

  G4double capBlackSheetZ[4] = {-WCBlackSheetThickness*zflip, 0., 0.,
                                                                (WCBarrelPMTOffset - (WCIDRadius-innerAnnulusRadius)) *zflip};
  G4double capBlackSheetRmin[4] = {0., 0., WCIDRadius, WCIDRadius};
  G4double capBlackSheetRmax[4] = {WCIDRadius+WCBlackSheetThickness,
                                   WCIDRadius+WCBlackSheetThickness,
                                                                   WCIDRadius+WCBlackSheetThickness,
                                                                   WCIDRadius+WCBlackSheetThickness};

Its placement within the plug volume is

  1828    G4VPhysicalVolume* physiWCCapBlackSheet =
  1829      new G4PVPlacement(0,
  1830                        G4ThreeVector(0.,0.,0.),
  1831                        logicWCCapBlackSheet,
  1832                        "WCCapBlackSheet",
  1833                        logicWCCap,
  1834                        false,
  1835                        0,
  1836                                            checkOverlaps);

z in plug

The placement at the origin within the WC[Top|Bot]Cap plug matches the black sheets internal coordinate system to that of the plug.

  z''_{BS,min} = -WCBlackSheetThickness
  z''_{BS,max} = WCBarrelPMTOffset - (WCIDRadius - innerAnnulusRadius)

In order for the blacksheet volume to be adjacent to border ring, it must have
z_{BS,max} == z_{plug,bev}
which does appear to be the case. Further, we must have
z_{BS,min} >= z_{plug,min}
which also seems to be the case.

After putting in some debugging output for the dimensions of the various geometry elements, i have found a problem. Here is an excerpt of that output that prints the dimensions computed for the WCTopCapBlackSheet:

WCTopCapBlackSheet dimensions:
 outer radius profile
    rmax[0] = 32420
    rmax[1] = 32420
    rmax[2] = 32420
    rmax[3] = 32420
 inner radius profile
    rmin[0] = 0
    rmin[1] = 0
    rmin[2] = 32400
    rmin[3] = 32400
 z profile (local coord)
    z'[0] = 20
    z'[1] = 0
    z'[2] = 0
    z'[3] = 1
 transform into TopCapAssembly (z of placement)
    z = z' + 2462.89
    z[0] = 2482.89
    z[1] = 2462.89
    z[2] = 2462.89
    z[3] = 2463.89

Note the z-profile in local coordinates:

 z profile (local coord)
    z'[0] = 20
    z'[1] = 0
    z'[2] = 0
    z'[3] = 1

z'[3] > z'[2] is a problem.

The values in the z profile vector are assumed to be monotonically increasing or decreasing (depending on zfilp). For the BlackSheet profile, this means that capBlackSheetZ[0] and capBlackSheetZ[3] should have opposite signs because capBlackSheetZ[1] == 0.. For the zflip = +1 case, capBlackSheetZ[0] < 0, so we should have

  capBlackSheetZ[3] >= 0
  WCBarrelPMTOffset - (WCIDRadius - innerAnnulusRadius) >= 0
  WCBarrelPMTOffset >= (WCIDRadius - innerAnnulusRadius)

This is a stronger restriction on WCBarrelPMTOffset than the condition for everything fitting into the CapAssembly volume, which was

  WCBarrelPMTOffset  >= WCIDRadius - innerAnnulusRadius - 1 mm

This explains why it was insufficient to add epsilon to WCBarrelPMTOffset. That was sufficient for the weaker condition, but insufficient for this opposite sign assumption made in the construction of the BlackSheet. Ok, so i need to add at least 1 mm to WCBarrelPMTOffset to satisfiy this strong condition.

Let me do that for the HyperK_HybridmPMT geometries and recheck. Oh, and while i am at it, put in a check for this condition. After this change, the z profile of the WCTopCapBlackSheet is now

 z profile (local coord)
    z'[0] = 20
    z'[1] = 0
    z'[2] = 0
    z'[3] = -0

Good!

Does that also resolve the overlaps? YES!

spradlin commented 1 year ago

OK, let us do another round of overlap checks. Running 1 event jobs of the following geometries with the WCSim_Geometry_Overlaps_CHECK build option turned ON:

Good. No more overlaps.

spradlin commented 1 year ago

I would like to bring this to a close by visualizing each of the endcaps. However, when i try to visualize the endcaps for any of the geometries that have an ExtraTower, the plug section of the cap does not display and i get the following error:

ERROR: G4VSceneHandler::RequestPrimitives
  Polyhedron not available for WCTopCap
  Touchable path: expHall 0 WCBox 0 WCBarrel 0 TopCapAssembly 0 WCTopCap 0 
  This means it cannot be visualized on most systems (try RayTracer).
  1) The solid may not have implemented the CreatePolyhedron method.
  2) For Boolean solids, the BooleanProcessor, which attempts to create
     the resultant polyhedron, may have failed.
ERROR: G4VSceneHandler::RequestPrimitives
  Polyhedron not available for WCTopCapBlackSheet
  Touchable path: expHall 0 WCBox 0 WCBarrel 0 TopCapAssembly 0 WCTopCap 0 WCTopCapBlackSheet 0 
  This means it cannot be visualized on most systems (try RayTracer).
  1) The solid may not have implemented the CreatePolyhedron method.
  2) For Boolean solids, the BooleanProcessor, which attempts to create
     the resultant polyhedron, may have failed.

For example, here is a visualization of the Top endcap of the HyperK_HybridmPMT_WithOD geomety without the ID PMTs displayed and with some of the barrel border cells set to invisible: HyperK_HybridmPMT_WithOD_TopCapAssembly_0001

I suspect that the Boolean construction of the cap plug and cap black sheet cannot be rendered. I do not know why. However, i think that should be followed up in another issue.

I also think that the investigation into the stopped particles should also be spun off into a separate issue.

And i think that the notes of this investigation into overlaps, and its accompanying merge request, is ready for review.

spradlin commented 1 year ago

This issue is closed with the merge of https://github.com/tdealtry/WCSim/pull/3.