[SOLVED]Trying to understand DistFromOutside

Hello,

I have a composite shape made of 4 rotated and translated half hyperboloids. This is suppose to model the intersection of 4 laser beams.

I’ve been working on a simulation of the interaction of a beam of atoms with these 4 laser beams. I do this using a montecarlo approach to decide if the atoms absorb photons from a particular beam (thus getting a momentum kick). They can then spontaneously decay (getting another kick). However, this simulation first has a free propagation part, which relies on the method DistFromOutside to check whether the atom will go through the volume defined by the beams and also estimate this distance.

Unfortunately, My program’s memory usages spikes incredibly when I shift the position of the laser beams causing some of the atoms to “miss” them (this does not happen when the atoms all go through the laser area, meaning that the method DistFromOutside gives a “finite” value).

I’m puzzled as to why this is happening. The following is a part of my code:

dist_atom_mol = beams->DistFromOutside( pos , dir );

if ( dist_atom_mol > 2 * molcham_rmin * 100 )

    atom_out = 1;

dtp[0] = ( dist_atom_mol / 100 ) / v.Mag() / sampling;


//------LVIS propagation prior to molasses beams------//

while ( atom_out == 0 && ! beams->Contains( pos ) )
{.....

The indicator atom_out just tells me if the atom will miss the laser region. Thus it will not enter the while loop that follows. I then use the computed distance to calculate a time step that I will use inside the loop.

I decided that I should better understand the DistFromOutside method to make sure that is not causing problems, so I wrote a very simple script (attached file) that creates a Box shape and then I create two points outside of it. Using DistFromOutside I ask what the distance to the box is. Point 1 correctly gives the distance. Point 2 distance to the box is intentionally asked in a direction that there is no possible way it will intersect the box. However, upon playing with the method’s options I don’t get the expected behavior.

For example if I choose iact = 0, I had no real expectations for what the value of safe_distance would be and I got 10. What bothered me is that despite specifying my own STEP = 100., the method still returned the default STEP value of 1.e30 (given by a call to some function Big() ).

Then I chose iact = 1 and I got the exact same behavior. I suppose Safe_distance = 10 (?). In line of this I lowered my proposed STEP to 9 so the method would return 9, as described in the manual. Well, it didn’t. It returned 1.e30 still and now safe_distance = 0. what!?

I tried options 2 and 3 and got same problems. I’ve tried all options on my actual code without any success. It’s possible I have some bugs lurking in my code and I will try to debug, but I was puzzled by DistFromOutside behavior.

Any help is greatly appreciated

Thanks!
create_box.C (1.63 KB)

Hi Andres,

The behavior of iact>0 in the computation of distances for shapes changed a bit due to consistency reasons for global navigation. Now any DistFromOutside that does not hit the shape return TGeoShape::Big() (meaning “no hit”) no matter the step proposed. Also, if the (partially) computed safety gets bigger than the proposed step, the final safety value is not computed for some shapes and a 0 value is returned.

This was just to explain the behavior of your simple macro. Now, you should note that using directly the navigation methods for solids with iact<3 (and even in general) is a bit for expert usage, this is why this is not described in the user manual. It is quite dangerous to use directly shape navigation methods since you can easily do mistakes. Imagine that you have a more complex hierarchic geometry and your interesting shape is somewhere in a branch, having a non-identical positioning matrix. Then all queries with respect to a point with global coordinates and direction have to be first converted to the local frame of the shape. Also, getting the distance and safety for a single shape may not be relevant if you have some other in the way of your track…

So, simply use global navigation queries from the level of TGeoManager, like: FindNode(), FindNextBoundary() or FindNextBoundaryAndStep(), Safety() or FindNormalFast(). These methods handle navigation in a global way, taking also into account rounding effects or boundary crossing issues.

Cheers,

Alright. Thanks for the advice. I’ve been trying to play with these to get the behavior I need but I just don’t know how to do it!

Before I was using DistFromOutside to get the distance each atom had to traverse to get to where the lasers are. I have now substituted this by FindNextBoundaryAndStep() + GetStep(). That is all good. However, I was also using the methods ‘Contains( pos )’ to tell me if the atom was or not in the laser region. Is using contains “safe”? or does it fall in the same category as DistFromOutside?

Furthermore, I need to know, in advance, to which volume does the boundary correspond. I guess this information is stored somehow in the pointer node that is output of FindNextBoundaryAndStep(), but how do I extract this information?

Thanks

Andres

After several trials I found a way using the GetPath() method (luckily I don’t have many nodes in my top volume). Though I had to use Step() and check the path after it, then retrace the Step. It’s not very elegant but seems to be working. I’ll do more tests.

Hi Andres,

Try FindNextBoundaryAndStep() which will also cross the next boundary. To find out where you are you can ask (besides GetPath()):
gGeoManager->GetCurrentNode(); // gets current node after crossing
gGeoManager-GetCurrentNode()->GetVolume(); // current volume

You can backup the volume pointer while building the geometry to check it during track propagation in order to avoid string comparison. Then the copy number of the node can give you which replica of the volume was hit.
Contains() at shape level is still not safe (see my previous post).

Taking the time to read the geometry chapter in the Users Guide may also save you some time :wink:
Cheers,

[quote=“agheata”]
Try FindNextBoundaryAndStep() which will also cross the next boundary.[/quote]

Yes this is what I’m using

Right now I’m doing string comparison as you said. This approach with volume seems good, I’ll give it a try.

[quote=“agheata”]
Taking the time to read the geometry chapter in the Users Guide may also save you some time [/quote]

I read it and then reread it, but I guess I’m just a newbie with ROOT and a not very experienced programmer.

Anyways, thanks a lot for your guidance.