rcs.f90

Fortran源文件

PROGRAM MAIN

integer,parameter::Imin=-40,Imax=40,Jmin=-40,Jmax=40,Kmin=-40,Kmax=40
integer,parameter::Iomin=Imin+5,Iomax=Imax-5,Jomin=Jmin+5,Jomax=Jmax-5,Komin=Kmin+5,Komax=Kmax-5
integer n
real*8 dx,z0,pi,lmd,wl,f,c,mu0,epsln0
real*8 k0,thetas,phis,RT1,RT2,rcs
complex*8 Jx,Jy,Jz,Jmx,Jmy,Jmz,CT1,CT2,ctemp

pi=3.1415926
c=3.e8
lmd=0.032
f=c/lmd
wl=42
dx=lmd/wl

epsln0=8.85e-12
mu0=pi*4.e-7
z0=sqrt(mu0/epsln0)
k0=2*pi/lmd

phis=90.
phis=phis*2*pi/180.
!thetas=0.
!thetas=thetas*2*pi/180.
Jx=cmplx(0,0)
Jy=cmplx(0,0)
Jz=cmplx(0,0)
Jmx=cmplx(0,0)
Jmy=cmplx(0,0)
Jmz=cmplx(0,0)

open(13,file='RCS1.dat')
open(14,file='RCS2.dat')
do n=0,180
  write(*,*)'Now the scattering angle thetas is:',n
  write(*,*)''
  thetas=n*2*pi/180.
  !phis=phis*2*pi/180.

  open(1,file='ExR.dat')
  open(2,file='ExI.dat')
  open(3,file='EyR.dat')
  open(4,file='EyI.dat')
  open(5,file='EzR.dat')
  open(6,file='EzI.dat')
  open(7,file='HxR.dat')
  open(8,file='HxI.dat')
  open(9,file='HyR.dat')
  open(10,file='HyI.dat')
  open(11,file='HzR.dat')
  open(12,file='HzI.dat')
  do j=Jomin,Jomax
    do k=Komin,Komax	  
	  ctemp=cmplx(0,k0*dx*((Iomin-0.5)*sin(thetas)*cos(phis)+j*sin(thetas)*sin(phis)+k*cos(thetas)))
      ctemp=cexp(ctemp)
	  read(3,*)RT1
	  read(4,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(5,*)RT1
	  read(6,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jmy=Jmy-CT2*ctemp
      Jmz=Jmz+CT1*ctemp
	  read(9,*)RT1
	  read(10,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(11,*)RT1
	  read(12,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jy=Jy+CT2*ctemp
      Jz=Jz-CT1*ctemp

	  ctemp=cmplx(0,k0*dx*((Iomax+0.5)*sin(thetas)*cos(phis)+j*sin(thetas)*sin(phis)+k*cos(thetas)))
      ctemp=cexp(ctemp)
	  read(3,*)RT1
	  read(4,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(5,*)RT1
	  read(6,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jmy=Jmy+CT2*ctemp
      Jmz=Jmz-CT1*ctemp
	  read(9,*)RT1
	  read(10,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(11,*)RT1
	  read(12,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jy=Jy-CT2*ctemp
      Jz=Jz+CT1*ctemp
	enddo
  enddo

  do i=Iomin,Iomax
    do k=Komin,Komax
	  ctemp=cmplx(0,k0*dx*(i*sin(thetas)*cos(phis)+(Jomin-0.5)*sin(thetas)*sin(phis)+k*cos(thetas)))
      ctemp=cexp(ctemp)
	  read(1,*)RT1
	  read(2,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(5,*)RT1
	  read(6,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jmx=Jmx+CT2*ctemp
      Jmz=Jmz-CT1*ctemp
	  read(7,*)RT1
	  read(8,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(11,*)RT1
	  read(12,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jx=Jx-CT2*ctemp
      Jz=Jz+CT1*ctemp
            
      ctemp=cmplx(0,k0*dx*(i*sin(thetas)*cos(phis)+(Jomax+0.5)*sin(thetas)*sin(phis)+k*cos(thetas)))
      ctemp=cexp(ctemp)
	  read(1,*)RT1
	  read(2,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(5,*)RT1
	  read(6,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jmx=Jmx-CT2*ctemp
      Jmz=Jmz+CT1*ctemp
	  read(7,*)RT1
	  read(8,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(11,*)RT1
	  read(12,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jx=Jx+CT2*ctemp
      Jz=Jz-CT1*ctemp
	enddo
  enddo

  do i=Iomin,Iomax
    do j=Jomin,Jomax
	  ctemp=cmplx(0,k0*dx*(i*sin(thetas)*cos(phis)+j*sin(thetas)*sin(phis)+(Komin-0.5)*cos(thetas)))
      ctemp=cexp(ctemp)
	  read(1,*)RT1
	  read(2,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(3,*)RT1
	  read(4,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jmx=Jmx-CT2*ctemp
      Jmy=Jmy+CT1*ctemp
	  read(7,*)RT1
	  read(8,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(9,*)RT1
	  read(10,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jx=Jx+CT2*ctemp
      Jy=Jy-CT1*ctemp
            
      ctemp=cmplx(0,k0*dx*(i*sin(thetas)*cos(phis)+j*sin(thetas)*sin(phis)+(Komax+0.5)*cos(thetas)))
      ctemp=cexp(ctemp)
	  read(1,*)RT1
	  read(2,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(3,*)RT1
	  read(4,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jmx=Jmx+CT2*ctemp
      Jmy=Jmy-CT1*ctemp
	  read(7,*)RT1
	  read(8,*)RT2
	  CT1=cmplx(RT1,RT2)
	  read(9,*)RT1
	  read(10,*)RT2
	  CT2=cmplx(RT1,RT2)
      Jx=Jx-CT2*ctemp
      Jy=Jy+CT1*ctemp
	enddo
  enddo
  close(1)
  close(2)
  close(3)
  close(4)
  close(5)
  close(6)
  close(7)
  close(8)
  close(9)
  close(10)
  close(11)
  close(12)

  Jmx=Jmx*dx*dx
  Jmy=Jmy*dx*dx
  Jmz=Jmz*dx*dx
  Jx=Jx*dx*dx
  Jy=Jy*dx*dx
  Jz=Jz*dx*dx
    
  rcs=abs(z0*(Jx*cos(thetas)*cos(phis)+Jy*cos(thetas)*sin(phis)-Jz*sin(thetas))+(-Jmx*sin(phis)+Jmy*cos(phis)))
  rcs=4*pi*(0.25*k0/pi*rcs)**2/lmd/lmd
  RCS=10*log10(rcs)
  write(13,*)n,RCS
  rcs=abs(z0*(Jx*sin(phis)-Jy*cos(phis))+(Jmx*cos(thetas)*cos(phis)+Jmy*cos(thetas)*sin(phis)-Jmz*sin(thetas)))
  rcs=4*pi*(0.25*k0/pi*rcs)**2/lmd/lmd
  RCS=10*log10(rcs)
  write(14,*)n,RCS
enddo
close(13)
close(14)

END