anchor.diff:  C----------------------------------------------------------------------C
   
         SUBROUTINE lstsqr(soln,gf,obs,nobs,OBSMAX,nflt,FLTMAX,damp,smooth,
  -     1                  smoof)
  +     1                  smoof,sgf,stscon,nsts)
   C----
   C Compute a least squares inversion for slip and rake
   C----
         IMPLICIT none
         CHARACTER*80 smoof
  +      REAL*8 stscon
  +      INTEGER nsts
         REAL*8 damp,smooth
         INTEGER OBSMAX,FLTMAX,dimx,dimy
         REAL*8 obs(OBSMAX,6)
  -      REAL*8 gf(OBSMAX,FLTMAX,6)
  -      REAL*8 A(3*OBSMAX+2*FLTMAX,2*FLTMAX)
  -      REAL*8 b(3*OBSMAX+2*FLTMAX,1)
  -      REAL*8 x(3*OBSMAX+2*FLTMAX,1)
  +      REAL*8 gf(OBSMAX,FLTMAX,6),sgf(FLTMAX,FLTMAX,4)
  +      REAL*8 A(3*OBSMAX+4*FLTMAX,2*FLTMAX)
  +      REAL*8 b(3*OBSMAX+4*FLTMAX,1)
  +      REAL*8 x(3*OBSMAX+4*FLTMAX,1)
         REAL*8 soln(FLTMAX,2)
         INTEGER i,j,nobs,nflt
         INTEGER n,m ! n=nrows, m=ncol
  
positive.code:               write(0,'(3F12.2,4F14.4)') (flt(i,j),j=1,7)
    902     continue
            write(0,*)
            write(0,'("HALF-SPACE PARAMETERS")')
            write(0,'("---------------------")')
            write(0,'(3A14)'),'vp','vs','dens'
            write(0,'(3F14.4)') haf(1),haf(2),haf(3)
            write(0,*)
        endif
  
  C----
  C Convert geographical to Cartesian coordinates, if necessary
  
negative.code: C Change number of rows, columns
        n = 3*nobs
        m = 2*nflt
  
  C----
  C Add stress constraints
  C----
        call addsts(A,b,n,m,sgf,stscon,nsts,OBSMAX,FLTMAX)
  
  C----
  C Add damping
  C----
        call damping(A,b,n,m,damp,OBSMAX,FLTMAX)
  
  C----
  C Add smoothing
  C----
        call smoothing(A,b,n,m,smooth,smoof,OBSMAX,FLTMAX)
  
  C----
  C Solve generalized least squares problem
  C----
        dimx = 3*OBSMAX+4*FLTMAX
        dimy = 2*FLTMAX
        call lsqsolve(A,b,x,dimx,dimy,n,m)
        do 404 i = 1,nflt
  