Print the parent particle information corresponding to the current particle

[chen2019 (chen)]

Dear Fluka experts,
I simulated a beam of photons hitting an object, and I wanted to print information about the current particle and the corresponding parent particle.I successfully printed the information of the current particle using USDRAW, but the information about the parent particle was not saved successfully. All the printed data was almost 0 (there were only two non-0 data in all the data, I don’t know why).

I consulted the manual and this post:
" One common application is the following: after an interaction which has produced sencondaries, let USDRAW copy some properties of the interacting particle into the TRACKR user variables. When STUPRF is called next to load the secondaries into stack, by default it copies the TRACKR user variables to the stack ones. In this way, information about the parent can be still carried by its daughters (and possibly by further descendants). This technique is sometimes referred to as “latching”."

Infomation about parent particle

Hello Antonino, In principle, you could save the information of the parent particles without using stupre.f and stuprf.f routines. When you have an interaction, entry USDRAW of mgdraw is called with the corresponding code of that type of interaction. Then, the properties of the secondaries produced in that interaction are stored in common GENSTK. At the moment of the interaction, you might want to save some of the properties of the parent in the TRACKR user variables SPAUSR(1-11) (for real v…

My understanding is: for the current particle, I print the information first and then save the information in TRACKR’s user variables. If the interaction produces secondary particles, STUPRF copies the user variables to the corresponding stack.When secondary particles start to transfer, the user variables in the stack are then copied into TRACKR’s user variables, and now I can get and print out information about the current particle and the parent particle of the previous cycle. Then I re-copy the information of the current particle as the parent particle to TRACKR.

      ENTRY USDRAW ( ICODE, MREG, XSCO, YSCO, ZSCO )
      IF ( .NOT. LFCOPE ) THEN
	  
         MRGNAM='DETUP'
         CALL GEON2R ( MRGNAM, NUPSTR, IERR )	  
         LFCOPE = .TRUE.

         WRITE ( 89, * ) '#Information about particles in region', MRGNAM
		 WRITE ( 89, * ) '#NCASE,ICODE,JTRACK,LTRACK, XSCO, YSCO, ZSCO,ETRACK'
		 WRITE ( 88, * ) '#Information about parent particles in region', MRGNAM
		 WRITE ( 88, * ) '#NCASE,ICODE,JTRACK,LTRACK, XSCO, YSCO, ZSCO,ETRACK'
      END IF

	  IF ( MREG .EQ. NUPSTR ) THEN
	     
*       current particles attributes		 
		 WRITE ( 89, * ) NCASE,ICODE,JTRACK,LTRACK, XSCO, YSCO, ZSCO,ETRACK
		 
*       the number of secondary particles and the corresponding particle ID		 
		 WRITE ( 89, * ) NP, (KPART(I),I=1,NP)
		 
*       parent particles attributes	  
		 WRITE ( 88, * ) ISPUSR(1),ISPUSR(2),ISPUSR(3),ISPUSR(4),SPAUSR(1),SPAUSR(2),SPAUSR(3),SPAUSR(4)
		 
*		Storing information about parent particles for secondary particles
		 ISPUSR(1) = NCASE
		 ISPUSR(2) = ICODE
         ISPUSR(3) = JTRACK
		 ISPUSR(4) = LTRACK
		 
		 SPAUSR(1) = XSCO
         SPAUSR(2) = YSCO
         SPAUSR(3) = ZSCO
		 SPAUSR(4) = ETRACK
		 	
	  END IF		 

      RETURN

mymgdraw_usdraw01.f (11.9 KB)

This is my USDRAW,I don’t know what causes it. Can someone help me analyze it？
Thanks

[dcalzola (Daniele Calzolari)]

Hello Chen,

your code is mostly right. Let me add a quick comment first. You correctly declare a logical initialization flag:

LOGICAL LFCOPE
DATA LFCOPE / .FALSE. /

Unfortunately, this flag is used both in the entries USDRAW and BXDRAW. Normally, you should enter first in the former, which then switches off the flag also for BXDRAW. This can explain why the .98 file is not appearing after the simulation.

Now let’s focus on the core of the problem. In the USDRAW entry you tell FLUKA to save the parent particle information in ISPUSR and SPAUSR variables only if the interaction took place in ‘DETUP’.
Therefore, all the interactions outside this region do not save any information.

To print out non-zero numbers in your .88 file, you would require two interactions in sequence in your geometry: the first one did not have any information saved (each variable is 0). After the first interaction, the secondary produced has the information of the parent particle. Now, if one of these secondaries interacts again in the same region, you save correctly the information.

To understand what happens in concrete with your simulation, I would need your inputfile to observe the geometry. My naive guess is that particles go into the detector area, interact once and then they leave.

I would suggest restructuring your routine saving the parent information outside of the if condition. Let me know if this solves your problems.

      IF ( MREG .EQ. NUPSTR ) THEN
     
*       current particles attributes		 
         WRITE ( 89, * ) NCASE,ICODE,JTRACK,LTRACK, XSCO, YSCO, ZSCO,ETRACK
		 
*       the number of secondary particles and the corresponding particle ID		 
         WRITE ( 89, * ) NP, (KPART(I),I=1,NP)
		 
*       parent particles attributes	  
         WRITE ( 88, * ) ISPUSR(1),ISPUSR(2),ISPUSR(3),ISPUSR(4),SPAUSR(1),SPAUSR(2),SPAUSR(3),SPAUSR(4)
		 
*		Storing information about parent particles for secondary particles
      END IF
      ISPUSR(1) = NCASE
      ISPUSR(2) = ICODE
      ISPUSR(3) = JTRACK
      ISPUSR(4) = LTRACK
      	 
      SPAUSR(1) = XSCO
      SPAUSR(2) = YSCO
      SPAUSR(3) = ZSCO
      SPAUSR(4) = ETRACK

As a personal final suggestion, I would encourage you to avoid using tabs and using just spaces. This, with a uniform indentation, strongly helps you in working with user routines.

Cheers,
Daniele

[chen2019 (chen)]

@dcalzola
Thank you for your comments.
I only use USDRAW in my routine because what(3) in my USERDUMP is set to 7.

@dcalzola

Now let’s focus on the core of the problem. In the USDRAW entry you tell FLUKA to save the parent particle information in ISPUSR and SPAUSR variables only if the interaction took place in ‘DETUP’.
Therefore, all the interactions outside this region do not save any information.

Yes, you’re right. That’s what I’m going to do.

@dcalzola

To print out non-zero numbers in your .88 file, you would require two interactions in sequence in your geometry: the first one did not have any information saved (each variable is 0). After the first interaction, the secondary produced has the information of the parent particle. Now, if one of these secondaries interacts again in the same region, you save correctly the information.

Yes, I have taken that into consideration and I can live with it. When primaries first interacts, the parent particle’s variables ISPUSR and SPAUSR do not hold values, so the print is 0. But the code should be able to record information about the parent particle and print it out as secondary particles continue to interact with second, third and more, and it shouldn’t be 0.

@dcalzola

To understand what happens in concrete with your simulation, I would need your inputfile to observe the geometry. My naive guess is that particles go into the detector area, interact once and then they leave.

As you can see from the LTRACK values in the fort.89 file, there are many values greater than 1, so there should still be a lot of secondary particle interactions occurring in the DETUP.
I don’t know if my comment is correct, this is my input file, can you see more problems.
Thanks
ScaSim111.flair (4.2 KB)
ScaSim111.inp (3.1 KB)

[dcalzola (Daniele Calzolari)]

As you can see from the LTRACK values in the fort.89 file, there are many values greater than 1, so there should still be a lot of secondary particle interactions occurring in the DETUP.

This is partially true. You can indeed see that some LTRACK are not 1, and this means that the particle is coming from a secondary interaction.

However, the photons do not interact a lot in the air of DETUP, and they start the interactions directly in the MUSCLE region. If a primary particle interacts there for the first time, all the secondaries will do not hold values in ISPUSR and SPAUSR, but they will still have LTRACK larger then one.

The very few values you observe are due to the photons interacting first in the DETUP region, and then going in the muscle region, finally going out again in the DETUP region and interacting another time there. (Or more simply, interacting twice in the DETUP region).

Now let’s go to the core problem of your simulation. You are using as a primary particles a 60 keV photons. By EMFCUT, you’re setting a transport cut of 1KeV for photons and 10 keV for electrons.
However, you did not specify nor the EMFCUT production cut nor the DELTARAY threshold.
I would suggest you to go as down as possible with the thresholds, selecting 100 eV for photons and 1 KeV for electrons.

What is happening then is the following: you have one photon being transported and interacting once.
All (or most of the produced particles) will fall below the production cut, and they are simply discarded.

Keep in mind that you are at the limit of the code, it might happen that very low energy electrons and photons are discarded before reaching the main loop.

[chen2019 (chen)]

@dcalzola
Thank you again for your patience and help.

@dcalzola

However, the photons do not interact a lot in the air of DETUP, and they start the interactions directly in the MUSCLE region. If a primary particle interacts there for the first time, all the secondaries will do not hold values in ISPUSR and SPAUSR, but they will still have LTRACK larger then one.

You are right, the DETUP region is filled air, the primary photon rarely interacts in this zone, and my code logically does not record the first interaction, so it prints all zeros.

So, as you suggested, I ran other tests. First, I added the production cutoff threshold, and the production and transmission cutoff thresholds for both photons and electrons were set to 1kev. Second, I tested it in the MUSCLE region,the printed data does change, generating a lot of data, but all the data is about ICODE=210.I checked the manual, ICODE=210 stands for Møller secondaries, but according to EEDRAW, does ICODE=210 still mean particles are below the production threshold?

Comparison of.89 and.88 files shows that, for example, primary particles 16 and 25, secondary particles from generation 1 to 8 have Compton and Rayleigh scattering in the MUSCLE area, but this information is not printed, only ICODE=210,why is that?
I’m sorry I don’t know much about this Møller secondaries. On the other hand, according to EEDRAW, this is probably due to the particle below production threshold.What does this have to do with my EMFCUT?
I apologize for the questions, but I would appreciate it if you have any comments
ScaSim222.inp (3.2 KB)
mymgdraw_usdraw01.f (12.0 KB)
ScaSim222001_fort.88.txt (1.1 MB)
ScaSim222001_fort.89.txt (1.2 MB)

[dcalzola (Daniele Calzolari)]

Do not apologize for your question, it is indeed extremely interesting and well-posed. Indeed, I checked the source code and it is not granted that in all cases the variables are copied from the TRACKR (ISPUSR, SPAUSR) to the EMFSTK. This can happen in some cases of electromagnetic interactions.

Without going into the technicalities of the code, since you are using very low-energy photons, the suggested approach would be the following:

1. You use stupre.f to save the parent particle values in EMFSTK.
2. As you are correctly doing, you continue to print these variables on a text file via mgdraw.f.

I attach here a possible stupre.f routine. Let me highlight some changes from the default behaviour:

1. I included the caslim.inc to have the information about NCASE
2. I manually selected the ICODE from the logic flags present in evtflg.inc file.

I attach the routine here. Try it with your case, and let me know if it works.
stupre.f (5.5 KB)

In case you are dealing with hadrons, muons or neutrinos, you would need to do a similar work with the stuprf.f routine. Let me know in case of further doubts.

Cheers, Daniele

[chen2019 (chen)]

Dear @dcalzola
I am very sorry for replying you so late,my PC was used to compute another simulation, which took many days, so I was late.
Thanks again for your help, the stupre.f routine code you provided seems to be working.The ISPUSR, SPAUSR variables I saved in usdraw were successfully printed in the .88 file.That’s what I’m trying to achieve, it’s perfect.
I learned your code and it seems that need to add code to stupre to manually save TRACKR information to EMFSTK. And, depending on your opinion, in my case, the default stupre routine may not do the job.
Anyway, thank you！