Photon backscattering depth measurement and modifying angular distribution graph

[jemancza (Jerzy Mikolaj Manczak)]

Hi Abdul,

Ok, now I think that we understand each other better. What you are trying to do is a good start, it would seem that you don’t need to modify the STUPRE function as the user flags should be propagated by default - in fact however, it has to be modified to be absolutely sure that the primary interaction data will be saved and propagated.
Please have a look at this thread:

Print the parent particle information corresponding to the current particle

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). [image] [image] I consulted the manual and this post: " One common applic…

Also, a discussion of a problem similar to yours can be found here:

Mgdraw routine to determine the number of Compton events before absorption

Hello David, Thank you for your feedback. I double checked the routine and indeed the problem was still related to flagging of some fluorescence photon created after Compton scattering of the primary photon. The (eventual) fluorescence photon and/or Auger electron were not correctly flagged due to the logic of the routine and in some rare cases the photon was inheriting the positive flag causing then the printing of its photoelectric effect from mgdraw. You can find here a stupre.f that now …

With your approach, you are asking for any Compton scattering, which would also have to be the very last interaction before your BXDRAW is invoked to see a non-zero SPAUSR(1).
The logic should like more like this:

ENTRY USDRAW ( ICODE, MREG, XSCO, YSCO, ZSCO )
IF (ICODE .EQ. 219 .AND. LTRACK.EQ.1) THEN
SPAUSR(1) = ZSCO
ENDIF

By setting a condition of LTRACK==1 you are ensuring that only the primary particle (in your case photons) interaction will modify the SPAUSR in the case where the very first interaction is a Compton scattering. Then you don’t touch the user variable anymore and the number you saved should be propagated until the very last generation. Now, you have to be careful here - what you are doing here will give you the ZSCO of the primary photons, but ONLY the ones whose first interaction was Compton scattering. What I was mentioning before in my previous post, you can in principle have a primary interaction which is different from Compton and still eventually get a signal in your detector.

THE ABOVE CODE WILL NOT WORK FOR YOUR PROBLEM, I just mentioned it to highlight what was conceptually wrong in the way you tried to approach it.

Due to the way that the FLUKA code is structured (you can find more details in the above forum discussions), you actually have to modify the stupre.f function and compile it with your own executable. Please find below two user routines that should do the job for you.
mgdraw_bdx_C2.f (7.0 KB)
stupre.f (4.1 KB)

Below I am discussing the most important code blocks of the attached files.
In the mgdraw_bdx_C2.f you have two non-empty functions:

      ENTRY BXDRAW ( ICODE, MREG, NEWREG, XSCO, YSCO, ZSCO )
	  IF ( .NOT. LFCOPE ) THEN
         LFCOPE = .TRUE.
         FILNAM1 = 'DetectedPhotons'       
         OPEN ( UNIT = 51, FILE = FILNAM1, STATUS = 'NEW')
      END IF

	  CALL GEOR2N ( MREG, MRGNAM, IERR1 )
	  CALL GEOR2N ( NEWREG, NRGNAM, IERR2 )
	  IF(IERR1 .NE. 0 .OR. IERR2 .NE. 0) STOP "Error in name conversion"
	  TFLAG = .FALSE.
	  IF(NRGNAM .EQ. "detector" .AND. MRGNAM.EQ."AIR") TFLAG = .TRUE. 
	  IF(TFLAG) WRITE(51,*) ICODE, NCASE, jtrack, ETRACK-AM(JTRACK), SPAUSR(1)
      RETURN

In this one, you save all the particles which cross the boundary between the air and the detector. It is important that NCASE is saved to associate the detected particles with the primary photons (more below). The SPAUSR(1) is the Z coordinate of the first interaction. Then you have:

      ENTRY USDRAW ( ICODE, MREG, XSCO, YSCO, ZSCO )   
	 IF ( .NOT. LFCOPE2 ) THEN
         LFCOPE2 = .TRUE.
         FILNAM2 = 'ScatteredPhotons'       
         OPEN ( UNIT = 52, FILE = FILNAM2, STATUS = 'NEW')
         END IF

         CALL GEOR2N ( MREG, MRGNAM, IERR1 )
	 IF(IERR1 .NE. 0 ) STOP "Error in name conversion"
         IF(MRGNAM .EQ. "TARGET" .AND. LTRACK.EQ.1) THEN
	    WRITE(52,*) ICODE, NCASE, jtrack, ETRACK-AM(JTRACK), XSCO, YSCO, ZSCO, MRGNAM
	 ENDIF

Where we save the information about the very first interaction of the primary photons (LTRACK == 1) whose interactions happened inside the target.
Now the magic happens in the stupre.f function:

      LOGICAL LBHABH
      LOGICAL PRIMER
      CHARACTER*8 MRGNAM
*
      CALL GEOR2N ( IREMF(1), MRGNAM, IERR1 )
      WRITE( LUNOUT,* ) MRGNAM
      IF(LTRACK.EQ.1 .AND. MRGNAM.EQ."TARGET") THEN
	 PRIMER = .TRUE.
      ELSE 
	 PRIMER = .FALSE.
      END IF

      IF ( LDLTRY ) THEN
         LBHABH = ICHEMF (NPEMF) * ICHEMF (NPEMF-1) .LT. 0
      ELSE
         LBHABH = .FALSE.
      END IF

      DO 500 NPNW = NPSTRT,NPEMF
	 IF (PRIMER) THEN 
	    ESPARK(1, NPNW) = ZEMF(NPNW)
	 ELSE
	    ESPARK(1, NPNW) = SPAUSR(1)
	 END IF
	 DO 51 ISPR = 2, MKBMX1
            ESPARK (ISPR,NPNW) = SPAUSR (ISPR)
 51      CONTINUE 
         DO 52 ISPR = 1, MKBMX2
            IESPAK (ISPR,NPNW) = ISPUSR (ISPR)
 52      CONTINUE
         LOUEMF (NPNW) = LLOUSE

The PRIMER flag is there to check whether we are dealing with the primary particle and the interaction happened inside the target. In fact, IREMF(1) is the region number of the first particle which is already loaded to the emf stack - stupre.f unlike stuprf.f does not have direct access to MREG - the last region of the track which is finishing its life.
The function action is different from the default one only when the PRIMER conditions are fulfilled - the last Z coordinate is saved to the user variable which then is correctly propagated until the end with the default behaviour of the STUPRE function.

Your simulation should finish with 2 files: ScatteredPhotons and DetectedPhotons. In ScatteredPhotons you have all the primary interactions with the information about the primary photon interaction. In DetectedPhotons you have the particles entering the detector together with the information of the primary interaction depth in Z. You can intersect both files with the NCASE variable.

ALTERNATIVELY

If you don’t want to mess with all the stupre.f magic, you can also skip this part and compile only the mgdraw_bdx_C2.f. Then you can still associate the primary interactions with the detected particles using NCASE and extract the information about the Z-depth only from ScatteredPhotons.

I hope that helps!

Cheers,
Jerzy