Skip to main content

Browse the archive

Revision 8cb2c8da847b4c6e53ee55172cf78036c95befd8 authored by M. Helena Gon\xe7alves on 06 June 2010, 00:00:00 UTC, committed by Gabor Csardi on 06 June 2010, 00:00:00 UTC 
None
0 parent
Raw File
mcpj.f 

      subroutine mcpj(th1,th2,psi1,p)
      implicit double precision(a-h,o-z)
      dimension p(2)

      ps1=psi1-1

C# calculate the probabilities of pj depending on psi1 pj=(p0,p1)

      if (dabs(ps1).gt.1.e-10) then

      d=dsqrt(1+ps1*(psi1*(th1-th2)**2-(th1+th2)**2+2*(th1+th2)))

        p(1)=(d-1+ps1*(th1-th2))/(2*ps1*(1-th2))

        p(2)=p(1)+(1-d+ps1*(th1+th2-2*th1*th2))/(2*ps1*th2*(1-th2))

      else
      p(1)=th1
      p(2)=th1

      end if
      return
      end
back to top