Production of charged pions by 730 MeV proton on thin lead target

Versions

FLUKA: 4.5.1
Flair: 3.4

Description

Dear colleagues,

I’d like to share with you the result for the comparison of the inclusive cross section for the production of charged pions by a 730 MeV/c proton on a thin (1-mm) lead target as predicted by FLUKA with experimental data. Experimental data is obtained from D. R. F. Cochran et al, Phys. Rev. D 6, 3085 (1972) Link

In the simulation, I score the production of all pions in the target by a customized mgdraw.f file, where I modified the USDRAW entry as follows:

   ENTRY USDRAW ( ICODE, MREG, XSCO, YSCO, ZSCO )
      IF (LFIRST) THEN
         LFIRST = .FALSE.
         IF ( KOMPUT .EQ. 2 ) THEN
            FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
         ELSE
            FILNAM = CFDRAW
         END IF
         PRINT *,'ENDRAW ',FILNAM
         CALL MYUSRINI(FILNAM)
c         OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
c     &          'FORMATTED' )
      END IF

c     Sample PRODUCTION of:
c     pim, pi0, pip 
c     photons
c     neutrino
      do ip=1,np
         if ((kpart(ip).eq.5).or.(kpart(ip).eq.6).or.
     &        (kpart(ip).eq.27).or.(kpart(ip).eq.28).or.
     &        (kpart(ip).eq.23).or.(kpart(ip).eq.13).or.
     &        (kpart(ip).eq.14).or.(kpart(ip).eq.7)) then
             CALL saveProd(ncase,kpart(ip),JTRACK,Tki(ip)+AM(kpart(ip)),
     &              Cxr(ip),Cyr(ip),Czr(ip),XSCO,YSCO,ZSCO,ETRACK)
          end if
       end do

* No output by default:
      RETURN

Here, the saveProd function is simply storing all the variables in a ROOT TTree(one entry for each call of saveProd).

To “compute” the FLUKA cross section, I consider, for each energy interval Delta E and cosine interval Delta c_z, the number of events therein Delta N. I am currently setting Delta E=10 MeV and Delta c_z=0.005.

The formula I am using to compute the cross section is:

dsigma/dE dOmega = (Delta N)(Delta E Delta cz ) * 1/(N0 * n * L) * 1/(2pi)

where L=0.1 cm is the target thickness, N0 is the total number of primaries in the simulation, and n=3.3E22/cm3 is the number density of lead atoms. The factor 2pi is to properly normalize the solid angle.

I’d like to share with you the results, to possibly get a feedback regarding the quality of the agreement.

I report the results in terms of the differential cross section as a function of the pion kinetic energy, for different values of the angle. An example is show in the following plot, for the pi+ production at an angle of 15 degrees. The data is reported as the white markers, while FLUKA results are in black / red. The black histogram is obtained from a fully analogue simulation (no bias at all), while the red is obtained from a simulation in which the interaction length of protons in the target is enhanced by a factor 10 via the LAM-BIAScard – in this case, the counts are scaled accordingly when computing the cross section.

All the results are in the attached PDF file.

out.pdf (197.1 KB)

Bests,

Andrea

Dear Andrea,
your FLUKA curves are correct.
(Note that it’s a 730 MeV kinetic energy proton beam, not 730 MeV/c momentum.)

Dear @ceruttif, thanks for the feedback. I apologize, it is a 730 kinetic-energy proton beam indeed :slight_smile: