32 ClassImp(TGeoShapeAssembly);
37 TGeoShapeAssembly::TGeoShapeAssembly()
47 TGeoShapeAssembly::TGeoShapeAssembly(TGeoVolumeAssembly *vol)
57 TGeoShapeAssembly::~TGeoShapeAssembly()
64 void TGeoShapeAssembly::ComputeBBox()
67 Fatal(
"ComputeBBox",
"Assembly shape %s without volume", GetName());
72 Int_t nd = fVolume->GetNdaughters();
73 if (!nd) {fBBoxOK = kTRUE;
return;}
76 Double_t xmin, xmax, ymin, ymax, zmin, zmax;
77 xmin = ymin = zmin = TGeoShape::Big();
78 xmax = ymax = zmax = -TGeoShape::Big();
81 for (Int_t i=0; i<nd; i++) {
82 node = fVolume->GetNode(i);
84 if (node->GetVolume()->IsAssembly()) node->GetVolume()->GetShape()->ComputeBBox();
85 box = (TGeoBBox*)node->GetVolume()->GetShape();
86 box->SetBoxPoints(vert);
87 for (Int_t ipt=0; ipt<8; ipt++) {
88 node->LocalToMaster(&vert[3*ipt], pt);
89 if (pt[0]<xmin) xmin=pt[0];
90 if (pt[0]>xmax) xmax=pt[0];
91 if (pt[1]<ymin) ymin=pt[1];
92 if (pt[1]>ymax) ymax=pt[1];
93 if (pt[2]<zmin) zmin=pt[2];
94 if (pt[2]>zmax) zmax=pt[2];
97 fDX = 0.5*(xmax-xmin);
98 fOrigin[0] = 0.5*(xmin+xmax);
99 fDY = 0.5*(ymax-ymin);
100 fOrigin[1] = 0.5*(ymin+ymax);
101 fDZ = 0.5*(zmax-zmin);
102 fOrigin[2] = 0.5*(zmin+zmax);
103 if (fDX>0 && fDY>0 && fDZ>0) fBBoxOK = kTRUE;
109 void TGeoShapeAssembly::RecomputeBoxLast()
111 Int_t nd = fVolume->GetNdaughters();
113 Warning(
"RecomputeBoxLast",
"No daughters for volume %s yet", fVolume->GetName());
116 TGeoNode *node = fVolume->GetNode(nd-1);
117 Double_t xmin, xmax, ymin, ymax, zmin, zmax;
119 xmin = ymin = zmin = TGeoShape::Big();
120 xmax = ymax = zmax = -TGeoShape::Big();
122 xmin = fOrigin[0]-fDX;
123 xmax = fOrigin[0]+fDX;
124 ymin = fOrigin[1]-fDY;
125 ymax = fOrigin[1]+fDY;
126 zmin = fOrigin[2]-fDZ;
127 zmax = fOrigin[2]+fDZ;
131 TGeoBBox *box = (TGeoBBox*)node->GetVolume()->GetShape();
132 if (TGeoShape::IsSameWithinTolerance(box->GetDX(), 0) ||
133 node->GetVolume()->IsAssembly()) node->GetVolume()->GetShape()->ComputeBBox();
134 box->SetBoxPoints(vert);
135 for (Int_t ipt=0; ipt<8; ipt++) {
136 node->LocalToMaster(&vert[3*ipt], pt);
137 if (pt[0]<xmin) xmin=pt[0];
138 if (pt[0]>xmax) xmax=pt[0];
139 if (pt[1]<ymin) ymin=pt[1];
140 if (pt[1]>ymax) ymax=pt[1];
141 if (pt[2]<zmin) zmin=pt[2];
142 if (pt[2]>zmax) zmax=pt[2];
144 fDX = 0.5*(xmax-xmin);
145 fOrigin[0] = 0.5*(xmin+xmax);
146 fDY = 0.5*(ymax-ymin);
147 fOrigin[1] = 0.5*(ymin+ymax);
148 fDZ = 0.5*(zmax-zmin);
149 fOrigin[2] = 0.5*(zmin+zmax);
156 void TGeoShapeAssembly::ComputeNormal(
const Double_t *point,
const Double_t *dir, Double_t *norm)
158 if (!fBBoxOK) ((TGeoShapeAssembly*)
this)->ComputeBBox();
159 Int_t inext = fVolume->GetNextNodeIndex();
161 DistFromOutside(point,dir,3);
162 inext = fVolume->GetNextNodeIndex();
164 Error(
"ComputeNormal",
"Invalid inext=%i (Ncomponents=%i)",inext,fVolume->GetNdaughters());
168 TGeoNode *node = fVolume->GetNode(inext);
169 Double_t local[3],ldir[3],lnorm[3];
170 node->MasterToLocal(point,local);
171 node->MasterToLocalVect(dir,ldir);
172 node->GetVolume()->GetShape()->ComputeNormal(local,ldir,lnorm);
173 node->LocalToMasterVect(lnorm,norm);
179 Bool_t TGeoShapeAssembly::Contains(
const Double_t *point)
const
181 if (!fBBoxOK) ((TGeoShapeAssembly*)
this)->ComputeBBox();
182 if (!TGeoBBox::Contains(point))
return kFALSE;
183 TGeoVoxelFinder *voxels = fVolume->GetVoxels();
186 Int_t *check_list = 0;
191 TGeoNavigator *nav = gGeoManager->GetCurrentNavigator();
192 TGeoStateInfo &td = *nav->GetCache()->GetInfo();
193 check_list = voxels->GetCheckList(point, ncheck, td);
195 nav->GetCache()->ReleaseInfo();
198 for (
id=0;
id<ncheck;
id++) {
199 node = fVolume->GetNode(check_list[
id]);
200 shape = node->GetVolume()->GetShape();
201 node->MasterToLocal(point,local);
202 if (shape->Contains(local)) {
203 fVolume->SetCurrentNodeIndex(check_list[
id]);
204 fVolume->SetNextNodeIndex(check_list[
id]);
205 nav->GetCache()->ReleaseInfo();
209 nav->GetCache()->ReleaseInfo();
212 Int_t nd = fVolume->GetNdaughters();
213 for (
id=0;
id<nd;
id++) {
214 node = fVolume->GetNode(
id);
215 shape = node->GetVolume()->GetShape();
216 node->MasterToLocal(point,local);
217 if (shape->Contains(local)) {
218 fVolume->SetCurrentNodeIndex(
id);
219 fVolume->SetNextNodeIndex(
id);
229 Int_t TGeoShapeAssembly::DistancetoPrimitive(Int_t , Int_t )
237 Double_t TGeoShapeAssembly::DistFromInside(
const Double_t * ,
const Double_t * , Int_t , Double_t , Double_t * )
const
239 Info(
"DistFromInside",
"Cannot compute distance from inside the assembly (but from a component)");
240 return TGeoShape::Big();
248 Double_t TGeoShapeAssembly::DistFromOutside(
const Double_t *point,
const Double_t *dir, Int_t iact, Double_t step, Double_t *safe)
const
253 TString sindent =
"";
254 for (Int_t k=0; k<indent; k++) sindent +=
" ";
255 Int_t idebug = TGeoManager::GetVerboseLevel();
257 if (!fBBoxOK) ((TGeoShapeAssembly*)
this)->ComputeBBox();
258 if (iact<3 && safe) {
259 *safe = Safety(point, kFALSE);
263 if (iact==0)
return TGeoShape::Big();
264 if ((iact==1) && (step<=*safe))
return TGeoShape::Big();
270 Double_t snext = 0.0;
272 Double_t stepmax = step;
275 Bool_t found = kFALSE;
276 memcpy(pt,point,3*
sizeof(Double_t));
278 if (idebug>4) printf(
"%s[%d] assembly %s checking distance to %d daughters...\n", sindent.Data(), indent, fVolume->GetName(), fVolume->GetNdaughters());
281 if (!TGeoBBox::Contains(point)) {
282 snext = TGeoBBox::DistFromOutside(point, dir, 3, stepmax);
284 snext = TMath::Min(0.01*snext, 1.E-6);
286 if (idebug>4 && snext > stepmax) printf(
"%s[%d] %s: bbox not crossed\n",sindent.Data(), indent, fVolume->GetName());
289 if (snext > stepmax)
return TGeoShape::Big();
293 for (i=0; i<3; i++) pt[i] += snext*dir[i];
306 Int_t nd = fVolume->GetNdaughters();
308 Double_t lpoint[3],ldir[3];
309 TGeoVoxelFinder *voxels = fVolume->GetVoxels();
310 if (nd<5 || !voxels) {
311 for (i=0; i<nd; i++) {
312 node = fVolume->GetNode(i);
313 if (voxels && voxels->IsSafeVoxel(pt, i, stepmax))
continue;
314 node->MasterToLocal(pt, lpoint);
315 node->MasterToLocalVect(dir, ldir);
317 if (idebug>4) printf(
"%s[%d] distance to %s ...\n", sindent.Data(), indent, node->GetName());
319 dist = node->GetVolume()->GetShape()->DistFromOutside(lpoint, ldir, 3, stepmax);
323 printf(
"%s[%d] %s -> from local=(%19.16f, %19.16f, %19.16f, %19.16f, %19.16f, %19.16f)\n",
324 sindent.Data(), indent, fVolume->GetName(), lpoint[0],lpoint[1],lpoint[2],ldir[0],ldir[1],ldir[2]);
325 printf(
"%s[%d] -> (l)to: %s shape %s snext=%g\n", sindent.Data(), indent, node->GetName(),
326 node->GetVolume()->GetShape()->ClassName(), dist);
331 fVolume->SetNextNodeIndex(i);
338 if (idebug>4) printf(
"%s[%d] %s: found %s at %f\n", sindent.Data(), indent, fVolume->GetName(), fVolume->GetNode(fVolume->GetNextNodeIndex())->GetName(), snext);
344 if (idebug>4) printf(
"%s[%d] %s: no daughter crossed\n", sindent.Data(), indent, fVolume->GetName());
347 return TGeoShape::Big();
352 TGeoNavigator *nav = gGeoManager->GetCurrentNavigator();
353 TGeoStateInfo &td = *nav->GetCache()->GetInfo();
355 voxels->SortCrossedVoxels(pt, dir, td);
356 while ((vlist=voxels->GetNextVoxel(pt, dir, ncheck, td))) {
357 for (i=0; i<ncheck; i++) {
358 node = fVolume->GetNode(vlist[i]);
359 node->MasterToLocal(pt, lpoint);
360 node->MasterToLocalVect(dir, ldir);
362 if (idebug>4) printf(
"%s[%d] distance to %s ...\n", sindent.Data(), indent, node->GetName());
364 dist = node->GetVolume()->GetShape()->DistFromOutside(lpoint, ldir, 3, stepmax);
368 printf(
"%s[%d] %s -> from local=(%19.16f, %19.16f, %19.16f, %19.16f, %19.16f, %19.16f)\n",
369 sindent.Data(), indent, fVolume->GetName(), lpoint[0],lpoint[1],lpoint[2], ldir[0],ldir[1],ldir[2]);
370 printf(
"%s[%d] -> to: %s shape %s snext=%g\n", sindent.Data(), indent, node->GetName(),
371 node->GetVolume()->GetShape()->ClassName(), dist);
375 fVolume->SetNextNodeIndex(vlist[i]);
380 nav->GetCache()->ReleaseInfo();
384 if (idebug>4) printf(
"%s[%d] %s: found %s at %f\n", sindent.Data(), indent, fVolume->GetName(), fVolume->GetNode(fVolume->GetNextNodeIndex())->GetName(), snext);
390 if (idebug>4) printf(
"%s[%d] %s: no daughter crossed\n", sindent.Data(), indent, fVolume->GetName());
393 return TGeoShape::Big();
399 TGeoVolume *TGeoShapeAssembly::Divide(TGeoVolume * ,
const char *divname, Int_t , Int_t ,
400 Double_t , Double_t )
402 Error(
"Divide",
"Assemblies cannot be divided. Division volume %s not created", divname);
410 TGeoShape *TGeoShapeAssembly::GetMakeRuntimeShape(TGeoShape * , TGeoMatrix * )
const
412 Error(
"GetMakeRuntimeShape",
"Assemblies cannot be parametrized.");
419 void TGeoShapeAssembly::InspectShape()
const
421 printf(
"*** Shape %s: TGeoShapeAssembly ***\n", GetName());
422 printf(
" Volume assembly %s with %i nodes\n", fVolume->GetName(), fVolume->GetNdaughters());
423 printf(
" Bounding box:\n");
424 if (!fBBoxOK) ((TGeoShapeAssembly*)
this)->ComputeBBox();
425 TGeoBBox::InspectShape();
431 void TGeoShapeAssembly::SetSegsAndPols(TBuffer3D & )
const
433 Error(
"SetSegsAndPols",
"Drawing functions should not be called for assemblies, but rather for their content");
440 Double_t TGeoShapeAssembly::Safety(
const Double_t *point, Bool_t in)
const
442 Double_t safety = TGeoShape::Big();
443 Double_t pt[3], loc[3];
444 if (!fBBoxOK) ((TGeoShapeAssembly*)
this)->ComputeBBox();
446 Int_t index = fVolume->GetCurrentNodeIndex();
447 TGeoVolume *vol = fVolume;
449 memcpy(loc, point, 3*
sizeof(Double_t));
451 memcpy(pt, loc, 3*
sizeof(Double_t));
452 node = vol->GetNode(index);
453 node->GetMatrix()->MasterToLocal(pt,loc);
454 vol = node->GetVolume();
455 index = vol->GetCurrentNodeIndex();
457 safety = vol->GetShape()->Safety(loc, in);
461 return TGeoShape::Big();
464 TGeoVoxelFinder *voxels = fVolume->GetVoxels();
465 Int_t nd = fVolume->GetNdaughters();
467 if (voxels) boxes = voxels->GetBoxes();
469 for (Int_t
id=0;
id<nd;
id++) {
473 Double_t dxyz0 = TMath::Abs(point[0]-boxes[ist+3])-boxes[ist];
474 if (dxyz0 > safety)
continue;
475 Double_t dxyz1 = TMath::Abs(point[1]-boxes[ist+4])-boxes[ist+1];
476 if (dxyz1 > safety)
continue;
477 Double_t dxyz2 = TMath::Abs(point[2]-boxes[ist+5])-boxes[ist+2];
478 if (dxyz2 > safety)
continue;
479 if (dxyz0>0) dxyz+=dxyz0*dxyz0;
480 if (dxyz1>0) dxyz+=dxyz1*dxyz1;
481 if (dxyz2>0) dxyz+=dxyz2*dxyz2;
482 if (dxyz >= safety*safety)
continue;
484 node = fVolume->GetNode(
id);
485 safe = node->Safety(point, kFALSE);
486 if (safe<=0.0)
return 0.0;
487 if (safe<safety) safety = safe;
495 void TGeoShapeAssembly::SavePrimitive(std::ostream & , Option_t * )
502 void TGeoShapeAssembly::SetPoints(Double_t * )
const
504 Error(
"SetPoints",
"Drawing functions should not be called for assemblies, but rather for their content");
510 void TGeoShapeAssembly::SetPoints(Float_t * )
const
512 Error(
"SetPoints",
"Drawing functions should not be called for assemblies, but rather for their content");
518 void TGeoShapeAssembly::GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols)
const
530 void TGeoShapeAssembly::Contains_v(
const Double_t *points, Bool_t *inside, Int_t vecsize)
const
532 for (Int_t i=0; i<vecsize; i++) inside[i] = Contains(&points[3*i]);
540 void TGeoShapeAssembly::ComputeNormal_v(
const Double_t *points,
const Double_t *dirs, Double_t *norms, Int_t vecsize)
542 for (Int_t i=0; i<vecsize; i++) ComputeNormal(&points[3*i], &dirs[3*i], &norms[3*i]);
548 void TGeoShapeAssembly::DistFromInside_v(
const Double_t *points,
const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step)
const
550 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromInside(&points[3*i], &dirs[3*i], 3, step[i]);
556 void TGeoShapeAssembly::DistFromOutside_v(
const Double_t *points,
const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t* step)
const
558 for (Int_t i=0; i<vecsize; i++) dists[i] = DistFromOutside(&points[3*i], &dirs[3*i], 3, step[i]);
566 void TGeoShapeAssembly::Safety_v(
const Double_t *points,
const Bool_t *inside, Double_t *safe, Int_t vecsize)
const
568 for (Int_t i=0; i<vecsize; i++) safe[i] = Safety(&points[3*i], inside[i]);