- 积分
- 13665
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2012-6-13
- 最后登录
- 1970-1-1
|
楼主 |
发表于 2013-5-4 10:16:27
|
显示全部楼层
c~~~~~~~~~~~~~~~~~~~~
parameter(m=252,n=97*37*2,kvt=2,n0=97*37)
dimension x(n,m),ect(kvt,m),v(n,kvt),oy0(n0,m),oy(n0,m),
$ ff0(n,m),ic(n0),x_ni(m,n)
open(2,file='e:\uv+\result\ecx.grd',form='binary',
$ access='direct',recl=n0*4)
irec1=1
do i=1,kvt
read(2,rec=irec1) (v(j,i),j=1,n)
irec1=irec1+1
write(*,*) v(7178,i)
enddo
irec=1
open(1,file='e:\uv-\u850_3.grd',
$form='binary',access='direct',recl=n0*4)
open(11,file='e:\uv-\v850_3.grd',
$form='binary',access='direct',recl=n0*4)
do ky=1,m
read(1,rec=irec) (oy0(i,ky),i=1,n0)
read(11,rec=irec) (oy(i,ky),i=1,n0)
irec=irec+1
do i=1,n0
if(abs(oy0(i,ky)).ge.999) then
ff0(i,ky)=-999.0
ic(i)=0
elseif(abs(oy(i,ky)).ge.999) then
ff0(i+n0,ky)=-999.0
ic(i)=0
else
ic(i)=1
ff0(i,ky)=oy0(i,ky)
ff0(i+n0,ky)=oy(i,ky)
endif
enddo
enddo
write(*,*) ff0(15,12)
call proj(n,m,kvt,v,ect,ff0)
open(111,file='e:\3\ect.txt')
do i=1,m
write(111,100) (ect(j,i),j=1,kvt)
100 format(1x,f10.4)
enddo
end
c~~~~~~~~~~~~~~~~~~~~~~~~~~~projection obs by V & eigenvalue~~~~~~~~~~~~~~
subroutine proj(nn,mm,kvt1,v1,ect1,x)
dimension x(nn,mm),ect1(kvt1,mm),v1(nn,kvt1),
$ v_ni(kvt1,nn)
do i=1,nn
do j=1,kvt1
v_ni(j,i)=v1(i,j)
enddo;enddo
do 50 it=1,mm
do 50 j=1,kvt1
ect1(j,it)=0.0
do 54 k=1,nn
54 ect1(j,it)=ect1(j,it)+v1(j,k)*x(k,it)
50 continue
END |
|