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

trapezoid rule was still wrong

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: 203
SVN UUID:     67feda4a-a26e-11df-9d6e-31afc202ad0c
parent bcc9f48d
c this is <greda.f>
c------------------------------------------------------------------------------
c $Id: greda.f,v 1.6 2000-06-21 21:02:21 thof Exp $
c $Id: greda.f,v 1.7 2000-06-23 21:47:20 thof Exp $
c
c 24/06/97 by Thomas Forbriger (IfG Stuttgart)
c
......@@ -66,6 +66,7 @@ c 24/05/00 V3.9 - changed cosine distance taper to apply factor greater
c than zero to last trace
c 21/06/00 V3.10 the trapezoid rule in subroutine backcoeff was awfully
c misbehaved
c 23/06/00 V3.11 trapezoid rule was still wrong by a constant factor 2
c
c==============================================================================
c
......@@ -74,7 +75,7 @@ c
c first we declare some general variables
c
character*79 version
parameter(version='GREDA V3.10 Greens function from data')
parameter(version='GREDA V3.11 Greens function from data')
c
c calculations common block
include 'greda.inc'
......@@ -1239,15 +1240,16 @@ c multiply with delta r to set up inverse transform
k=j
j=chain(k)
if (l.eq.k) then
dr=0.5*(r(j)-r(k))
dr=(r(j)-r(k))
elseif (j.lt.1) then
dr=0.5*(r(k)-r(l))
dr=(r(k)-r(l))
else
dr=(r(j)-r(l))
endif
dr=dr*0.5
alpha(k)=alpha(k)*dr
if (verbose) print *,' offset ',i,r(k),'m step ',dr,'m',
& ' array index ',k
if (verbose) print *,' offset ',i,r(k),'m trapezoid step ',dr,'m',
& ' index ',k
enddo
c
c print numbers
......
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