Message: Re: Primative Scorers in Parameterised Volume Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Primative Scorers in Parameterised Volume 

Keywords: Primative Scorer, Hits, Voxel Phantom
Forum: Hits, Digitization and Pileup
Re: Question Primative Scorers in Parameterised Volume (Jonathan)
Re: None Re: Primative Scorers in Parameterised Volume (Makoto Asai)
Date: 05 Jul, 2006
From: Jonathan <Jonathan>

Hi Makoto,

Thanks for your response. Below is my DetectorConstruction, followed by the parameterisation where I place the voxels:

DetectorConstruction:

void VoxelPhantomGeometry::PhantomConstructor()
{
	VoxelPhantomImporter* importPhantom = new VoxelPhantomImporter;
	Voxel voxels;
	G4int cols,rows,slices = 0;
	G4double xDim,yDim,zDim = 0;

	G4String fileName = "";	

	G4cout << "Please enter the name of the segmented voxel file (enter 1 for default): " << G4endl;
	G4cin >> fileName;

	if(fileName == "1")
        fileName = "examplePatient.csv";

	importPhantom->ReadASCIIFile(fileName,voxels);

	cols = voxels.getVoxelCols();
	rows = voxels.getVoxelRows();
	slices = voxels.getVoxelSlices();

	G4ThreeVector dim = voxels.getVoxelDimensions();
	xDim = dim.x();
	yDim = dim.y();
	zDim = dim.z();

    //Building up the parameterisation ...

    G4Box* voxelParameterisation = new G4Box( "Voxels", (xDim/2)*mm, 
    							(yDim/2)*mm, 
    							(zDim/2)*mm);

    voxelParameterisationLogicalvolume = new G4LogicalVolume(voxelParameterisation,
    					                     phantom,
    						             "VoxelParameterisation_logical");

    G4int numberOfVoxels = voxels.getTotalNumberOfVoxels();

    G4VPVParameterisation* paramVoxel = new VoxelPhantomParameterisation(lunginhale,
						                         lungexhale,
				                                         adiposeTissue,
						                         breast,
							                 phantom,
									 muscle,
								         liver,
									 denseBone,
								         trabecularBone );

    voxelParameterisationPhysical = new G4PVParameterised( "Voxels (physical)", 
                            					voxelParameterisationLogicalvolume, 
                            					physiWorld,
                            					kUndefined, 	// 3D optimization
                            					numberOfVoxels,  // copy's
                            					paramVoxel );			
    delete importPhantom;

}

void VoxelPhantomGeometry::SetupDetectors()
{
    //================================================
    // Sensitive detectors : MultiFunctionalDetector
    //================================================

    G4SDManager* SDman = G4SDManager::GetSDMpointer();

    G4String voxelSDname = "voxelSD";

    //------------------------
    // MultiFunctionalDetector
    //------------------------

    G4MultiFunctionalDetector* MFDet = new G4MultiFunctionalDetector(voxelSDname);
    SDman->AddNewDetector( MFDet );                 // Register SD to SDManager.
    voxelParameterisationLogicalvolume->SetSensitiveDetector(MFDet);  // Assign SD to the logical volume.

    //---------------------------------------
    // SDFilter : Sensitive Detector Filters
    //---------------------------------------

    G4String fltName,particleName;

    //-- gamma filter
    G4SDParticleFilter* gammaFilter = 
        new G4SDParticleFilter(fltName="gammaFilter", particleName="gamma");

    //-- electron filter
    G4SDParticleFilter* electronFilter = 
        new G4SDParticleFilter(fltName="electronFilter");
    electronFilter->add(particleName="e+");   // accept electrons.
    electronFilter->add(particleName="e-");   // accept positorons.

    //------------------------
    // PS : Primitive Scorers
    //------------------------

    G4String psName;

    // total
    G4PSEnergyDeposit*      scorer0 = new G4PSEnergyDeposit(psName="totalEDep", 0);
    G4PSDoseDeposit*        scorer1 = new G4PSDoseDeposit(psName="totalDose", 0);

    // by gammas
    G4PSEnergyDeposit*      scorer2 = new G4PSEnergyDeposit(psName="gammaEDep", 0);
    scorer2->SetFilter(gammaFilter);

    G4PSDoseDeposit*        scorer3 = new G4PSDoseDeposit(psName="gammaDose", 0);
    scorer3->SetFilter(gammaFilter);

    // by electrons    
    G4PSEnergyDeposit*      scorer4 = new G4PSEnergyDeposit(psName="electronEDep", 0);
    scorer4->SetFilter(electronFilter);

    G4PSDoseDeposit*        scorer5 = new G4PSDoseDeposit(psName="electronDose", 0);
    scorer5->SetFilter(electronFilter);

    //------------------------------------------------------------
    //  Register primitive scorers to MultiFunctionalDetector
    //------------------------------------------------------------

    MFDet->RegisterPrimitive(scorer0);
    MFDet->RegisterPrimitive(scorer1);
    MFDet->RegisterPrimitive(scorer2);
    MFDet->RegisterPrimitive(scorer3);
    MFDet->RegisterPrimitive(scorer4);
    MFDet->RegisterPrimitive(scorer5);	
}

G4VPhysicalVolume* VoxelPhantomGeometry::Construct()
{
	InitialisationOfMaterials();

    G4double worldXDimension = 5000.*mm;
    G4double worldYDimension = 5000.*mm;
    G4double worldZDimension = 5000.*mm;

    solidWorld = new G4Box( "WorldSolid",
                            worldXDimension,
                            worldYDimension,
                            worldZDimension );

    logicWorld = new G4LogicalVolume( solidWorld, 
                                      air, 
                                      "WorldLogical", 
                                      0, 0, 0 );

    physiWorld = new G4PVPlacement( 0,
                                    G4ThreeVector(0,0,0),
                                    "World",
                                    logicWorld,
                                    0,
                                    false,
                                    0 );

    //logicWorld->SetVisAttributes (G4VisAttributes::Invisible);	   

    PhantomConstructor();

    SetupDetectors();

    return physiWorld;	

}

And my Parameterisation: (xPos, yPos, zPos declared static)

void VoxelPhantomParameterisation::ComputeTransformation(const G4int copyNo, 
													G4VPhysicalVolume* physVol) const
{	
      if(copyNo == 0)
      {
	xPos = -colNumber/2.;
	yPos = -rowNumber/2.;
	zPos = -sliceNumber/2.;
      } 

 	  // increments row
      if((copyNo % colNumber == 0) && copyNo != 0)
      {
      	xPos = -colNumber/2.;
      	yPos++;
      }

      // increments slice
      if((copyNo % (colNumber * rowNumber) == 0) && copyNo != 0) 
	  {
	  	xPos = -colNumber/2.;
	  	yPos = -rowNumber/2.;
	  	zPos++;
	  }

	  // increments column
	  if((copyNo % colNumber != 0) && (copyNo % (colNumber * rowNumber) !=0) && copyNo != 0)
	  {
	  	xPos++;
	  } 

	  // positions all voxels so that their origin will fall in line with mother's origin (0,0,0)	  
	  G4ThreeVector origin(((xPos*xDim)+(xDim/2.))*mm, 
	  					   ((yPos*yDim)+(yDim/2.))*mm, 
	  					   ((zPos*zDim)+(zDim/2.))*mm);
          // debuggin
	  //G4cout << origin << ", ";

	  physVol->SetTranslation(origin);

	  // prep for an additional run (ie, if the copyNo is equal to the last element, reset all counters
	  if((copyNo + 1) == (colNumber * rowNumber * sliceNumber))
	  {
	  	xPos = -colNumber/2.;
	  	yPos = -rowNumber/2.;
	  	zPos = -sliceNumber/2.;
	  }

}

void VoxelPhantomParameterisation::ComputeDimensions(G4Box& voxels, const G4int, 
													 const G4VPhysicalVolume*) const
{
      //not needed yet until you want to update on the fly or change geometry as function of copyNo
      voxels.SetXHalfLength(voxelDim.x()/2 * mm);
      voxels.SetYHalfLength(voxelDim.y()/2 * mm);
      voxels.SetZHalfLength(voxelDim.z()/2 * mm);
} 

G4Material* VoxelPhantomParameterisation::ComputeMaterial(const G4int copyNo, G4VPhysicalVolume* physVol,
	                                                	  const G4VTouchable*)
{
      // uses the material map formed in VoxelImporter class to map the materials to the voxels
	  G4int chooser = materials[copyNo];

	  switch(chooser)
	  {
	  	case 1:          physVol->SetName("PhysicalLungINhale");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeLungINhale );
                         return lungInhale;
	  	case 11:          physVol->SetName("PhysicalAdipose");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeAdipose );
                         return adiposeTissue;
	  	case 0:          physVol->SetName("PhysicalBreast");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeBreast );
                         return breastTissue;
	  	case 4:          physVol->SetName("PhysicalPhantom");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributePhantom );
                         return phantomTissue;
	  	case 14:          physVol->SetName("PhysicalMuscle");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeMuscle );
                         return muscleTissue;
	  	case 35:          physVol->SetName("PhysicalLiver");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeLiver );
                         return liverTissue;
	  	case 10:          physVol->SetName("PhysicalTrabecularBone");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeTrabecularBone );
                         return trabecularBoneTissue;
	  	case 8:          physVol->SetName("PhysicalDenseBone");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributeDenseBone );
                         return denseBoneTissue;                         
	  	default: 	  	 physVol->SetName("PhysicalPhantom");
                         physVol->GetLogicalVolume()->SetVisAttributes( attributePhantom );
                         return phantomTissue;
	  }

return physVol->GetLogicalVolume()->GetMaterial(); }

Thanks again for taking the time to help me out on this. There is also a related post that I made that may shed some light on the problem in the Events and Tracking forum #509. If you need anymore info, just let me know!

Jonathan

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

1 None: Re: Primative Scorers in Parameterised Volume   (Sylvia Studeny - 05 Jul, 2006)
1 None: Re: Primative Scorers in Parameterised Volume   (Jonathan - 06 Jul, 2006)
3 None: Re: Primative Scorers in Parameterised Volume   (Jonathan - 06 Jul, 2006)
 Add Message Add Message
to: "Re: Primative Scorers in Parameterised Volume"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]