barnes滤波

yeah...

利用 barnes进行空间滤波 fortran程序一直运算但是生成不了数据  求大神解答

program main
implicit none
integer,parameter::xx=480,yy=241,zz=27,tt=48,dg=1
real   ,parameter :: undef=-32767
!global scale,1*1
real,parameter::C1=4000,C2=16000,G=0.3,pi=3.14

real,external::distance
integer::h,i,j,k,l,m,n
real::r
real::hgt(xx,yy,zz,tt),F1(xx,yy,zz,tt),F2(xx,yy,zz,tt),F(xx,yy,zz,tt),D1(xx,yy,zz,tt),D2(xx,yy,zz,tt),E1(xx,yy,zz,tt),E2(xx,yy,zz,tt),F11(xx,yy,zz,tt),F21(xx,yy,zz,tt),F12(xx,yy,zz,tt),F22(xx,yy,zz,tt),t1,t2,w1,w2

CHARACTER*70 :: FILEPATH1='G:\lunwen_0.75degree\0601\data\ec\decomposional\ud0601_0.75degree\'
CHARACTER*130 :: FILEPATH2='G:\lunwen_0.75degree\0601\data\ec\decomposional\smooth\barnes_smooth\ud0601_0.75degree_c1-4000c2-16000\'

!!!!!!!!!!!!!!!!!!!!!!!!!!!    read initial field        !!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! open(10,file='E:\ScientificResearchDuringMaster\Wrf_out\19-21\wavefiltering\0719PRESwavefilter.dat',access='direct',recl=4,status='old')
! open(20,file='E:\ScientificResearchDuringMaster\Wrf_out\19-21\wavefiltering\PRESfiltered.grd',form='binary')
! open(30,file='E:\ScientificResearchDuringMaster\Wrf_out\19-21\wavefiltering\PRESditong.grd',form='binary')
OPEN(10,FILE=trim(FILEPATH1)//'\bin\outud0601_0.75degree20060508_20060519_quyu.bin',FORM='binary',ACCESS='DIRECT',recl=4,status='old')
OPEN(20,FILE=trim(FILEPATH2)//'\bin\outud0601_0.75degree20060508_20060519_quyu_baners.bin',FORM='binary')
OPEN(30,FILE=trim(FILEPATH2)//'\bin\outud0601_0.75degree20060508_20060519_quyu_baner_ditong.bin',FORM='binary')

  h=0
do m=1,tt
print *,'Begin time = ',m
do n=1,zz
  do j=1,yy
   do i=1,xx
    h=h+1
    read(10,rec=h) hgt(i,j,n,m)
    print *, hgt(1,10,20,40)
enddo
enddo
enddo
enddo

do m=1,tt
do n=1,zz
  do j=1,yy
   do i=1,xx
    w1=0
    t1=0
    w2=0
    t2=0
   do k=1,yy,dg
    do l=1,xx,dg
     r=distance(-90.0+(j-1)*0.75,-90.0+(k-1)*0.75,(i-1)*0.75,(l-1)*0.75)
     if(r<=sqrt(C1))then
      w1=w1+exp(-r**2/(4*C1))
      t1=t1+hgt(l,k,n,m)*exp(-r**2/(4*C1))
      print*,r
     endif
     if(r<=sqrt(C2))then
      w2=w2+exp(-r**2/(4*C2))
      t2=t2+hgt(l,k,n,m)*exp(-r**2/(4*C2))
     endif
enddo
enddo
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
   F1(i,j,n,m)=t1/w1     !meso scale
   F2(i,j,n,m)=t2/w2     !large scale
enddo
enddo
enddo
enddo

D1=hgt-F1
D2=hgt-F2

do m=1,tt
do n=1,zz
  do j=1,yy
   do i=1,xx
   w1=0
   t1=0
   w2=0
   t2=0
   do k=1,yy,dg
    do l=1,xx,dg
     r=distance(-90.0+(j-1)*0.75,-90.0+(k-1)*0.75,(i-1)*0.75,(l-1)*0.75)
     if(r<=sqrt(C1))then
      w1=w1+exp(-r**2/(4*G*C1))
      t1=t1+D1(l,k,n,m)*exp(-r**2/(4*G*C1))
     endif
     if(r<=sqrt(C2))then
      w2=w2+exp(-r**2/(4*G*C2))
      t2=t2+D2(l,k,n,m)*exp(-r**2/(4*G*C2))
     endif
    enddo
   enddo

   F11(i,j,n,m)=F1(i,j,n,m)+t1/w1
   F21(i,j,n,m)=F2(i,j,n,m)+t2/w2
enddo
enddo
enddo
enddo
   E1=F11-F1
   E2=F21-F2

do m=1,tt
do n=1,zz
  do j=1,yy
   do i=1,xx
   w1=0
   t1=0
   w2=0
   t2=0
   do k=1,yy,dg
    do l=1,xx,dg
     r=distance(-90.0+(j-1)*0.75,-90.0+(k-1)*0.75,(i-1)*0.75,(l-1)*0.75)
     if(r<=sqrt(C1))then
      w1=w1+exp(-r**2/(4*C1))
      t1=t1+E1(l,k,n,m)*exp(-r**2/(4*C1))
     endif
     if(r<=sqrt(C2))then
      w2=w2+exp(-r**2/(4*C2))
      t2=t2+E2(l,k,n,m)*exp(-r**2/(4*C2))
     endif
    enddo
   enddo
   F12(i,j,n,m)=F11(i,j,n,m)+3*(F11(i,j,n,m)-F1(i,j,n,m))/4-t1/w1
   F22(i,j,n,m)=F21(i,j,n,m)+3*(F21(i,j,n,m)-F2(i,j,n,m))/4-t2/w2
   F(i,j,n,m)=F12(i,j,n,m)-F22(i,j,n,m)
   write(20) F(i,j,n,m)
   write(30) F22(i,j,n,m)
   print *, i,j,n,m
enddo
enddo
enddo

print *,'End time = ',m
enddo


close(10)
close(20)
close(30)
end program main

! write(*,*) hgt(1:10,1:10)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!   calculate the 'F0' in the Barnes formula!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!  calculate the 'F2' under different 'C' in the Barnes formula (the second revision)  !!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!calculate the spherical distance between (i,j)  and (k,l)!!!!!!!!!!!!!!!!
function distance(lat1,lat2,lon1,lon2)
implicit none
real,parameter::R=6371.004,pi=3.1415926
real::lat1,lat2,lon1,lon2
real::distance
distance=R*acos(cos(lat1*pi/180)*cos(lat2*pi/180)*cos((lon1-lon2)*pi/180)+sin(lat1*pi/180)*sin(lat2*pi/180))
return
end function distance

ljchen1989

二维滤波不是一般用loess么？

yeah...

ljchen1989 发表于 2017-2-15 20:47
二维滤波不是一般用loess么？

您说的是用matlab的loess进行空间滤波吗  主要是不太会matlab

ljchen1989

yeah... 发表于 2017-2-15 21:43
您说的是用matlab的loess进行空间滤波吗  主要是不太会matlab

会matlab也不一定会loess滤波啊~
那神器，会了直接发nature

李白

非常感谢，赞！！！

夏夜

ljchen1989 发表于 2017-2-15 20:47
二维滤波不是一般用loess么？

loess二维滤波怎么实现？有么有资料啊

pop_lxh

解决了吗，我也想要程序

jiuzhou526

谢谢  大赞