Export volumes of given material from TGeoManager

Dear Experts,
I am using TGeoManager::Import(“aaa.gdml”) to get a geometry description. I would like to keep only the volumes that contain specific material and export the geometry.
It seems possible to loop from the top volume and check what material is in each sub volume but I couldn’t find a way to remove a volume.
Do you see any solution?
Thank in advance for your help.

Maybe @agheata can give some hints

Hi,
Just keeping the volumes wih some material is not enough because they are connected as a tree via their TGeoNode daughter list, so you may need to keep also volumes not having the same materials bu having daughters that do have.
You can try using the TGeoIterator crawler (see for example tutorials/geom/runplugin.C) and use TGeoVolume::RemoveNode to remove the geometry branches you are not interested in. However, you may need to do 2 passes with the crawler to mark all volumes to be kept first. You may have some problems because the geometry as imported from gdml is closed (operation fixing stuff and building optimization structures) and cannot be re-opened, so you will need to call manually some methods like TGeoVolume::Voxelize to have the chopped geometry valid. Exporting again may work though. Good luck!

2 Likes

Thanks for the idea @agheata
I fact, I have first exported the GDML to a .C macro so I can insert the filtering before the geometry is closed. It is working with a simple .C file.
With the .C file created by gGeoManager->Export("bla.C"), I have lots of errors when I try to read the file with root bla.C. There are many redefinition of variables or improper syntax like:

pMed1 = new TGeoMedium(...)

instead of

 TGeoMedium* pMed1 = new TGeoMedium(...)

Do you know how to fix this?

I have managed to filter the nodes according to their material but when I export the geometry, all the node are still there. I think it is because the geometry is already closed. @agheata do you see a way to reopen it?

I think I found a solution. I basically copied the sequence to import the gdml (see ROOT: TGeoManager Class Reference) file but removing the closing of the geometry that I do after the filtering.

{

  // import from a gdml file without closing the geometry
  TGeoManager* GeoManager = new TGeoManager("GDMLImport", "Geometry imported from GDML");
  TString cmd = TString::Format("TGDMLParse::StartGDML(\"%s\")", "NA62.gdml");
  TGeoVolume* world = (TGeoVolume*)gROOT->ProcessLineFast(cmd);  
  GeoManager->SetTopVolume(world);


  /// cleaning
  TFile *_file0 = TFile::Open("NA62Nu.root","RECREATE");
  TGeoIterator next( GeoManager->GetTopVolume () );

  TString path;
  TGeoNode *node;
  //  vector<TString> keeps;
  set<TString> keep;
  while((node=next())) {
    next.GetPath(path);
    if(TString(node->GetVolume()->GetMaterial()->GetName()) == TString("G4_lKr") ||
       TString(node->GetVolume()->GetMaterial()->GetName()) == TString("G4_Fe")
       ){
      	cout<<"Keep "<<path<<endl;
      TObjArray *tx = path.Tokenize("/");
      for (Int_t i = 0; i < tx->GetEntries(); i++){
	//	keeps.push_back( ((TObjString *)(tx->At(i)))->String() );
	keep.insert( ((TObjString *)(tx->At(i)))->String() );
      }
    }
  }
  
  while( 1 ){
    next.Reset(GeoManager->GetTopVolume ());
    while((node=next())) {
      TString name = node->GetName();
      if ( keep.find(name) == keep.end()){
	node->GetMotherVolume()->RemoveNode(node);
	cout<<"Remove "<<name<<endl;
	next.Reset(GeoManager->GetTopVolume ());
	break;
      }
    }
    if( !(node=next()) ){
      cout<<"Nothing else to remove"<<endl;
      break;
    }
  }

  next.Reset(GeoManager->GetTopVolume ());
  while((node=next())) {
    

  }
  GeoManager->CloseGeometry();
  GeoManager->Write();
  _file0->Close();
}





I still have a problem. While the RemoveNode() command is indeed removing the node from the geometry, the volumes on the contrary are not removed.
@agheata is there a way to do this?

There is no support for this and volumes may be cross-referenced, but you can try to use (after removing all unnecessary nodes):

   auto volumes = gGeoManager->GetListOfVolumes();
   // loop volumes to remove
   volumes->Remove(one_volume);
   // end loop
   volumes->Compress();
2 Likes

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.