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

debugging was successful

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: 2288
SVN UUID:     67feda4a-a26e-11df-9d6e-31afc202ad0c
parent 08dcf62a
c this is <gresynoise.f>
c ----------------------------------------------------------------------------
c ($Id: gresynoise.f,v 1.9 2007-05-22 11:51:59 tforb Exp $)
c ($Id: gresynoise.f,v 1.10 2007-05-24 14:38:21 tforb Exp $)
c
c Copyright (c) 2007 by Thomas Forbriger (BFO Schiltach)
c
......@@ -13,9 +13,8 @@ c 09/05/2007 V1.0 Thomas Forbriger
c 10/05/2007 V1.1 introduced randomization
c 11/05/2007 V1.2 allow selection of number of output samples
c 15/05/2007 V1.3 reordered loops; phase is of correct type now
c 21/05/2007 V1.4i provide longer seismograms by keep multiple copies of
c 21/05/2007 V1.4 provide longer seismograms by keep multiple copies of
c noisymized Fourier coefficients
c not yet finished
c
c ============================================================================
c
......@@ -46,10 +45,10 @@ c
c
character*(*) version
parameter(version=
&'GRESYNOISE V1.4i calculate noise seismograms')
&'GRESYNOISE V1.4 calculate noise seismograms')
character*(*) GRESYNOISE_CVS_ID
parameter(GRESYNOISE_CVS_ID=
&'$Id: gresynoise.f,v 1.9 2007-05-22 11:51:59 tforb Exp $')
&'$Id: gresynoise.f,v 1.10 2007-05-24 14:38:21 tforb Exp $')
c
c dimensions
integer maxtr, maxsetsamp, maxom, maxu, maxset, maxsamp
......@@ -115,19 +114,20 @@ c options
logical debug, optlambda, optnew, optresponse, hankel1, hankel2
logical suppress, optrndscale
real lambdalim, tapfrac, scalexp
integer limnsetsamp
integer argnsets
c commandline
integer maxopt, lastarg, iargc
parameter(maxopt=11)
parameter(maxopt=12)
character*2 optid(maxopt)
character*80 optarg(maxopt)
logical optset(maxopt), opthasarg(maxopt)
c here are the keys to our commandline options
data optid/2h-d,2h-l,2h-o,2h-t,2h-i,2h-S,2h-1,2h-2,2h-e,
& 2h-r,2h-n/
& 2h-r,2h-n,2h-m/
data opthasarg/.FALSE.,.TRUE.,.FALSE.,.TRUE.,4*.FALSE.,
& .TRUE.,.FALSE.,.TRUE./
data optarg/1h-,2h1.,1h-,3h10.,4*1h-,2h0.,2*1h1/
& .TRUE.,.FALSE.,2*.TRUE./
data optarg/1h-,2h1.,1h-,3h10.,4*1h-,2h0.,3*1h1/
c------------------------------------------------------------------------------
c basic information
c
......@@ -136,7 +136,7 @@ c
if ((greenname(1:5).eq.'-help').or.(iargc().lt.5)) then
print *,version
print *,'Usage: gresynoise Zgreen Rgreen Zseis Rseis rcvfile'
print *,' [-l lambda] [-o] [-n n]'
print *,' [-l lambda] [-o] [-n n] [-m n]'
print *,' [-t frac] [-e exp] [-r]'
print *,' [-S] [-1] [-2] [-i]'
print *,'or gresynoise -help'
......@@ -174,6 +174,8 @@ c
print *,'-n n number of sets to use'
print *,' time series will become larger with each'
print *,' additional set of Fourier coefficients'
print *,'-m n limit number of samples in one set to n'
print *,' if possible'
print *,' '
print *,'Input file units:'
print *,' '
......@@ -214,6 +216,8 @@ c
argnsets=max(1,argnsets)
print *,'number of sets to be used: ',argnsets
nset=argnsets
limnsetsamp=maxsetsamp
if (optset(12)) read(optarg(12), *, err=99) limnsetsamp
c
pi2=2.d0*pi
c
......@@ -226,6 +230,16 @@ c
print *,'Rseisname: ',Rseisname(1:index(Rseisname, ' ')-1)
print *,'rcvname: ',rcvname(1:index(rcvname, ' ')-1)
endif
if (debug) then
print *,'array dimensions:'
print *,'maxu: ',maxu
print *,'maxtr: ',maxtr
print *,'maxom: ',maxom
print *,'maxset: ',maxset
print *,'maxsamp: ',maxsamp
print *,'maxsetsamp: ',maxsetsamp
endif
c
c----------------------------------------------------------------------
c read configuration and greens function
......@@ -337,11 +351,14 @@ c
& stop 'ERROR: I assume the first frequency to be 0.Hz!'
c we have got all frequencies now and are going to fit them to
c a power of 2 number
limnsetsamp=limnsetsamp/2
limnsetsamp=min(maxsetsamp/2,limnsetsamp)
limnsetsamp=max(Znom,limnsetsamp)
pown=0
1 continue
pown=pown+1
nsetsamp=2**pown
if (nsetsamp.lt.maxsetsamp) goto 1
if (nsetsamp.lt.limnsetsamp) goto 1
if (nsetsamp.lt.(2*Znom)) stop 'ERROR: cannot use enough samples'
c calculate sampling interval
dt=pi2/(nsetsamp*dom)
......@@ -452,7 +469,7 @@ c scale trace and stack
scal=scal*scalefac(i,is)
cc print *,'rec # ', i, ' scal: ',scal
do io=1,Znom
print *,i,io,is,Zfourier(io),Rfourier(io),scal
cc print *,i,io,is,Zfourier(io),Rfourier(io),scal
Zsets(io,is)=Zsets(io,is)+Zfourier(io)*scal
Rsets(io,is)=Zsets(io,is)+Rfourier(io)*scal
enddo
......@@ -473,7 +490,7 @@ c
c
do is=1,nset
c make Fourier coefficients complete
do io=1,maxom
do io=1,nsetsamp
Ztrans(io)=(0.d0,0.d0)
Rtrans(io)=(0.d0,0.d0)
enddo
......@@ -489,9 +506,17 @@ c
print *,'transform set ',is,' to time domain...'
c transform to time domain
if (debug) print *,'go fork'
cc print *,'go fork'
cc do io=1,nsetsamp
cc print *, io, Ztrans(io), Rtrans(io)
cc enddo
call tf_dfork(nsetsamp, Ztrans, 1.d0)
call tf_dfork(nsetsamp, Rtrans, 1.d0)
if (debug) print *,'done fork'
cc print *,'done fork'
cc do io=1,nsetsamp
cc print *, io, Ztrans(io), Rtrans(io)
cc enddo
c extract time series
ltap=nsetsamp/2
rtap=ltap+1
......
Supports Markdown
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