PROGRAM Compute_SZProfiles ! A code for computing: ! 1st column: x=r/r500 ! 2nd column: SZ profile in RJ temperature [uK] ! using Vikhlinin's model parameters. ! E. Komatsu, November 17, 2011 IMPLICIT none double precision :: xin=0d0,xout=6d0 double precision :: y,Tsz,tol=1d-6 double precision :: n0,rc,rs,alpha,beta,gamma,eps,n02,rc2,beta2,r500v,t3,et3 double precision :: para(13),pe,pe2d character(len=10) :: name integer :: i,icls,ncls=76 double precision :: rombint external pe open(2,file='name.txt',status='old') do icls=1,ncls read(2,*)name open(3,file='vikhlinin_xray_parameters/'//name,status='old') read(3,*)n0,rc,rs,alpha,beta,gamma,eps,n02,rc2,beta2,r500v,t3,et3 print*,name print*,'r500=',r500v*0.72d0/1d3,' h^-1 Mpc' print*,'pressure profile is cut at r/r500=',xout open(1,file='x_tsz_'//trim(name)//'.txt') do i=1,1000 y=(xin+1d-3*dble(i)*(xout-xin)) ! y=theta/theta500 para=(/y,n0,rc,rs,alpha,beta,gamma,eps,n02,rc2,beta2,r500v,t3/) pe2d=(r500v/1d3)*rombint(pe,-dsqrt(xout**2d0-y**2d0),dsqrt(xout**2d0-y**2d0),tol,13,para) Tsz=-2d0*273d0*(pe2d/25d0) ! uK write(1,'(2E18.5)')y,tsz enddo close(1) enddo close(2) END PROGRAM Compute_SZProfiles !----------------------------------------------------------------------------- DOUBLE PRECISION FUNCTION pe(l,y) ! x=r/r500 and r500 is in units of kpc ! pe = electron pressure IMPLICIT none double precision :: l,x,y(13) double precision :: n0,rc,rs,alpha,beta,gamma,eps,n02,rc2,beta2,r500,t3 double precision :: r,nesq,te x=dsqrt(l**2d0+y(1)**2d0) n0=y(2) rc=y(3) rs=y(4) alpha=y(5) beta=y(6) gamma=y(7) eps=y(8) n02=y(9) rc2=y(10) beta2=y(11) r500=y(12) t3=y(13) Te=1.35d0*(t3/1.11d0/0.97d0)*((x/0.045d0)**1.9d0+0.45d0)/((x/0.045d0)**1.9d0+1d0) & /(1d0+(x/0.6d0)**2d0)**0.45d0 r=x*r500 ! kpc nesq=n0**2d0*(r/rc)**(-alpha)/(1d0+r**2d0/rc**2d0)**(3d0*beta-alpha/2d0) & /(1d0+r**gamma/rs**gamma)**(eps/gamma) + n02**2d0/(1d0+r**2d0/rc2**2d0)**(3d0*beta2) pe=dsqrt(nesq)*(te*1d3) ! eV cm^-3, electron pressure return END FUNCTION pe