zhangyiss/DSurfTomo
1
0
Fork 0
mirror of https://github.com/HongjianFang/DSurfTomo.git synced 2025-05-08 00:01:14 +08:00
Code Issues Actions Packages Projects Releases Wiki Activity

runable for group velocity

Browse Source 
still need to be tested for Love wave group vel
This commit is contained in:
Hongjian Fang 2017-07-03 05:51:18 +08:00
parent db87c3c9d8
commit f4181ebfa7
4 changed files with 19 additions and 17 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

11
scripts/GenerateIniMOD.py
View File

@ -6,19 +6,20 @@ import numpy as np
#parameters need to be changed #parameters need to be changed
#start #start
nx=18 nx=75
ny=18 ny=96
nz=8 nz=17
minvel=0.8 minvel=0.8
velgrad=0.5 velgrad=0.5
dep1=np.array([0,0.2,0.4,0.6,0.8,1.1,1.4,1.8,2.5]) dep1=1.5+np.array([-1.5, 0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0,10.0,11.0,13.0,16.0,20.0,30.0])
vel1=np.loadtxt('mod.1d')
#end #end
vs1=np.zeros(nz) vs1=np.zeros(nz)
mod=np.zeros((nz*ny,nx)) mod=np.zeros((nz*ny,nx))
for k in range(nz): for k in range(nz):
for j in range(ny): for j in range(ny):
for i in range(nx): for i in range(nx):
mod[k*ny+j,i]= minvel+dep1[k]*velgrad mod[k*ny+j,i]= vel1[k]/1.75#minvel+dep1[k]*velgrad
with open('MOD','w') as fp: with open('MOD','w') as fp:
for i in range(nz): for i in range(nz):
fp.write('%5.1f' % dep1[i]) fp.write('%5.1f' % dep1[i])

23
src/CalSurfG.f90
View File

@ -973,7 +973,7 @@ subroutine CalSurfG(nx,ny,nz,nparpi,vels,iw,rw,col,dsurf, &
real vpz(nz),vsz(nz),rhoz(nz),depz(nz) real vpz(nz),vsz(nz),rhoz(nz),depz(nz)
real*8 pvRc(nx*ny,kmaxRc),pvRg(nx*ny,kmaxRg),pvLc(nx*ny,kmaxLc),pvLg(nx*ny,kmaxLg) real*8 pvRc(nx*ny,kmax),pvRg(nx*ny,kmaxRg),pvLc(nx*ny,kmax),pvLg(nx*ny,kmaxLg)
real*8 sen_vsRc(nx*ny,kmaxRc,nz),sen_vpRc(nx*ny,kmaxRc,nz) real*8 sen_vsRc(nx*ny,kmaxRc,nz),sen_vpRc(nx*ny,kmaxRc,nz)
real*8 sen_rhoRc(nx*ny,kmaxRc,nz) real*8 sen_rhoRc(nx*ny,kmaxRc,nz)
real*8 sen_vsRg(nx*ny,kmaxRg,nz),sen_vpRg(nx*ny,kmaxRg,nz) real*8 sen_vsRg(nx*ny,kmaxRg,nz),sen_vpRg(nx*ny,kmaxRg,nz)
@ -1076,8 +1076,9 @@ subroutine CalSurfG(nx,ny,nz,nparpi,vels,iw,rw,col,dsurf, &
if(kmaxRg.gt.0) then if(kmaxRg.gt.0) then
iwave=2 iwave=2
igr=0 igr=0
! print*,kmax
call caldespersion(nx,ny,nz,vels,pvRc, & call caldespersion(nx,ny,nz,vels,pvRc, &
iwave,igr,kmaxRc,tRc,depz,minthk) iwave,igr,kmax,tRg,depz,minthk)
igr=1 igr=1
call depthkernel(nx,ny,nz,vels,pvRg,sen_vsRg,sen_vpRg, & call depthkernel(nx,ny,nz,vels,pvRg,sen_vsRg,sen_vpRg, &
sen_rhoRg,iwave,igr,kmaxRg,tRg,depz,minthk) sen_rhoRg,iwave,igr,kmaxRg,tRg,depz,minthk)
@ -1087,14 +1088,14 @@ subroutine CalSurfG(nx,ny,nz,nparpi,vels,iw,rw,col,dsurf, &
iwave=1 iwave=1
igr=0 igr=0
call depthkernel(nx,ny,nz,vels,pvLc,sen_vsLc,sen_vpLc, & call depthkernel(nx,ny,nz,vels,pvLc,sen_vsLc,sen_vpLc, &
sen_rhoLc,iwave,igr,kmaxLc,tLc,depz,minthk) sen_rhoLc,iwave,igr,kmax,tLc,depz,minthk)
endif endif
if(kmaxLg.gt.0) then if(kmaxLg.gt.0) then
iwave=1 iwave=1
igr=0 igr=0
call caldespersion(nx,ny,nz,vels,pvLc, & call caldespersion(nx,ny,nz,vels,pvLc, &
iwave,igr,kmaxRc,tRc,depz,minthk) iwave,igr,kmax,tLg,depz,minthk)
igr=1 igr=1
call depthkernel(nx,ny,nz,vels,pvLg,sen_vsLg,sen_vpLg, & call depthkernel(nx,ny,nz,vels,pvLg,sen_vsLg,sen_vpLg, &
sen_rhoLg,iwave,igr,kmaxLg,tLg,depz,minthk) sen_rhoLg,iwave,igr,kmaxLg,tLg,depz,minthk)
@ -1372,11 +1373,17 @@ subroutine CalSurfG(nx,ny,nz,nparpi,vels,iw,rw,col,dsurf, &
enddo enddo
endif ! 'if' before rpath endif ! 'if' before rpath
enddo enddo
IF(asgr.EQ.1)THEN
DEALLOCATE (velnb, STAT=checkstat)
IF(checkstat > 0)THEN
WRITE(6,*)'Error with DEALLOCATE: PROGRAM fmmin2d: velnb'
ENDIF
ENDIF
IF(asgr.EQ.1)DEALLOCATE(ttnr,nstsr)
enddo ! 'do' before gridder enddo ! 'do' before gridder
IF(asgr.EQ.1)DEALLOCATE(ttnr,nstsr)
IF(rbint.EQ.1)THEN IF(rbint.EQ.1)THEN
WRITE(6,*)'Note that at least one two-point ray path' WRITE(6,*)'Note that at least one two-point ray path'
@ -1386,12 +1393,6 @@ subroutine CalSurfG(nx,ny,nz,nparpi,vels,iw,rw,col,dsurf, &
WRITE(6,*)'that you adjust the dimensions of your grid' WRITE(6,*)'that you adjust the dimensions of your grid'
WRITE(6,*)'to prevent this from occurring.' WRITE(6,*)'to prevent this from occurring.'
ENDIF ENDIF
IF(asgr.EQ.1)THEN
DEALLOCATE (velnb, STAT=checkstat)
IF(checkstat > 0)THEN
WRITE(6,*)'Error with DEALLOCATE: PROGRAM fmmin2d: velnb'
ENDIF
ENDIF
enddo enddo
enddo enddo
deallocate(fdm) deallocate(fdm)

BIN
src/DSurfTomo
View File

Binary file not shown.

2
src/main.f90
View File

@ -498,7 +498,7 @@ program SurfTomo
acond = 0.0 acond = 0.0
arnorm = 0.0 arnorm = 0.0
xnorm = 0.0 xnorm = 0.0
localSize = n/4 localSize = 10
call LSMR(m, n, leniw, lenrw,iw,rw,cbst, damp,& call LSMR(m, n, leniw, lenrw,iw,rw,cbst, damp,&
atol, btol, conlim, itnlim, localSize, nout,& atol, btol, conlim, itnlim, localSize, nout,&