anchor.diff:  C   1 for finite source
   C----
         call readfaults(nflt,evlo,evla,evdp,str,dip,rak,slip,dy,dx,area,
  -     1                flttyp,maxflt)
  +     1                flttyp,FMAX)
   C----
   C Calculate net displacements at each station
   C----
  
positive.code: C an isotropic elastic half-space due to shear dislocations.
  C----
        IMPLICIT NONE
        LOGICAL verbose
        CHARACTER*20 srcfile,stafile,haffile
        INTEGER nflt,maxflt
        PARAMETER (maxflt=1500)
        INTEGER flttyp(maxflt)
        REAL*8 evlo(maxflt),evla(maxflt),evdp(maxflt),str(maxflt),
       1       dip(maxflt),rak(maxflt),slip(maxflt),dx(maxflt),dy(maxflt),
       2       area(maxflt)
  
        REAL*8 stlo,stla,stdp
        REAL*8 vp,vs,dens
  
        REAL*8 uNnet,uEnet,uZnet
  
        call gcmdln(verbose)
  C----
  C Input files
  C----
        call checkctrlfiles(srcfile,stafile,haffile,verbose)
        open (unit=21,file=srcfile,status='old')
        open (unit=22,file=stafile,status='old')
        open (unit=23,file=haffile,status='old')
        read (23,*) vp,vs,dens
  C----
  C Output file
  C----
  
negative.code:           goto 101
    103 continue
  
   9999 format (2(f10.3),3(f14.4))
  
        END
  
  C======================================================================C
  
        SUBROUTINE readfaults(nflt,evlo,evla,evdp,str,dip,rak,slip,dy,dx,
       1                      area,flttyp,FMAX)
  
        IMPLICIT none
        INTEGER flt,nflt,FMAX,flttyp(FMAX)
        REAL*8 evlo(FMAX),evla(FMAX),evdp(FMAX),str(FMAX),
       1       dip(FMAX),rak(FMAX),slip(FMAX),dy(FMAX),dx(FMAX),
       2       area(FMAX)
  
        flt = 0
    201 flt = flt + 1
        read (21,*,END=202) evlo(flt),evla(flt),evdp(flt),str(flt),
  