Commit 281c7a91 authored by thomas.forbriger's avatar thomas.forbriger Committed by thomas.forbriger
Browse files

introduced line source

This is a legacy commit from before 2015-03-01.
It may be incomplete as well as inconsistent.
See COPYING.legacy and README.history for details.


SVN Path:     http://gpitrsvn.gpi.uni-karlsruhe.de/repos/TFSoftware/trunk
SVN Revision: 2997
SVN UUID:     67feda4a-a26e-11df-9d6e-31afc202ad0c
parent df39d45b
c this is <refmet.f> originally by J. Ungerer 1990
c======================================================================
c $Id: refmet.f,v 1.12 2009-02-23 18:22:29 tforb Exp $
c $Id: refmet.f,v 1.13 2010-03-02 13:29:18 tforb Exp $
c
c Reflectivity Method
c
......@@ -101,14 +101,15 @@ c 29/04/2000 V2.7 introduced Hankel functions
c 11/02/2003 V2.8 introduced new parameter definition (for const-Q case
c only): tabluated velocity refers to real part of modulus
c 23/02/2009 V2.9 porting code to gfortran
c 02/03/2010 V1.10 implementing line source
c
c======================================================================
PROGRAM refmet
character*70 version
parameter(version='REFMET V2.9 Reflectivity Method')
parameter(version='REFMET V2.10 Reflectivity Method')
character*79 cvsid
parameter(cvsid='$Id: refmet.f,v 1.12 2009-02-23 18:22:29 tforb Exp $')
parameter(cvsid='$Id: refmet.f,v 1.13 2010-03-02 13:29:18 tforb Exp $')
c array dimension declaration
integer me, msl, mf, ms
......@@ -618,6 +619,7 @@ c 29/04/00 * introduced H0, H1, H2 for Hankel function evaluation
c * and the switches cl_hankel1 and cl_hankel2
c * declared Neumann functions
c * introduced cl_hankelswitch
c 02/03/10 * cl_linesrc
c
c======================================================================
c used libraries and external subroutines
......@@ -728,10 +730,11 @@ c verbosity
c Hankel functions
double complex H0, H1, H2
logical cl_hankel1, cl_hankel2
logical cl_linesrc
double precision cl_hankelswitch
c commandline
integer cl_maxopt, cl_lastarg, iargc
parameter(cl_maxopt=8)
parameter(cl_maxopt=9)
character*2 cl_optid(cl_maxopt)
character*40 cl_optarg(cl_maxopt)
logical cl_optset(cl_maxopt), cl_opthasarg(cl_maxopt)
......@@ -740,9 +743,9 @@ c command line options and parameter values
character*80 cl_comsel
logical cl_debug, cl_rprogress, cl_comeach
c here are the keys to our commandline options
data cl_optid/2h-d,2h-v,2h-o,2h-c,2h-s,2h-1,2h-2,2h-p/
data cl_opthasarg/.FALSE.,2*.TRUE.,.FALSE.,.TRUE.,2*.FALSE.,.TRUE./
data cl_optarg/1h-,1h1,10hrefmet.out,1h-,2h ,2*1h-,2h0./
data cl_optid/2h-d,2h-v,2h-o,2h-c,2h-s,2h-1,2h-2,2h-p,2h-l/
data cl_opthasarg/.FALSE.,2*.TRUE.,.FALSE.,.TRUE.,2*.FALSE.,.TRUE.,.FALSE./
data cl_optarg/1h-,1h1,10hrefmet.out,1h-,2h ,2*1h-,2h0.,1h-/
C**********************************************************************C
C 1. Herdfunktion (als Funktionsdefinition) C
......@@ -799,6 +802,7 @@ c set options
cl_hankel1=cl_optset(6)
cl_hankel2=cl_optset(7)
read(cl_optarg(8), *) cl_hankelswitch
cl_linesrc=cl_optset(9)
cl_rprogress=.false.
if (cl_vlevel.lt.0) cl_rprogress=.true.
cl_vlevel=abs(cl_vlevel)
......@@ -1305,6 +1309,8 @@ c slowness range
sff_free(sff_nfree)='use the first Hankel functions H^(1)_m'
elseif (cl_hankel2) then
sff_free(sff_nfree)='use the second Hankel functions H^(2)_m'
elseif (cl_linesrc) then
sff_free(sff_nfree)='use line source (cosine kernel)'
else
sff_free(sff_nfree)='use the Bessel functions J_m'
endif
......@@ -2091,6 +2097,37 @@ c no tangential component, no near field
IDNz=dcmplx(0.d0)
endif
elseif (cl_linesrc)
C Line source calculation with cosine kernel
J0=cos(jarg)
J1=0.
J2=0.
c calculate integrands
if (typ.eq.1) then
stop 'ERROR: line source not defined for moment tensor'
elseif (typ.eq.2) then
c for a vertical single force
IDFr=dcmplx(0.d0)
c only vertical component is defined for line source
if (u.gt.1.e-50) then
IDFz=2.d0*FFI*IME*In4*J0*w(f)/u
else
IDFz=dcmplx(0.d0)
endif
c no tangential component, no near field
IDFphi=dcmplx(0.d0)
IDNphi=dcmplx(0.d0)
IDNr=dcmplx(0.d0)
IDNz=dcmplx(0.d0)
endif
endif
else
C Besselfunktionen aufrufen
J0=d_j0(jarg)
......
......@@ -73,6 +73,12 @@ c
print *,' Bessel function J_m for wavefield expansion.'
print *,'-p u use Hankel functions only for slowness values'
print *,' greater than ''u'' (in s/km).'
print *,'-l use line source (cosine kernels)'
print *,' note: only vertical single force is defined'
print *,' radial component is not defined and will'
print *,' be set to zero.'
print *,' For technical reasons the integration will'
print *,' not include slowness euqal zero.'
print *,' '
print *,'file Is the name of the main configuration file.'
print *,' It contains the names of the three file'
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment