anchor.diff:  !----
   ! Input fault depth and dimensions depend on GF model
   if (fault%file.ne.'none') then
  +
  +    ! Rectangular faults in an elastic half-space (okada_rect)
       if (gf_model.eq.'okada_rect') then
  +
  +        ! Depth is positive down, units in meters
           if (minval(fault%array(:,3)).lt.0.0d0) then
  -            ! Depth is positive down
  +            ! Found negative depth
               write(stderr,*) 'read_inputs: found fault depth of ',minval(fault%array(:,3))
  -            call usage('Depth is defined positive down for gf_model '//trim(gf_model))
  +            call usage('Depth is defined positive down for gf_model '//trim(gf_model)// &
  +                       '(usage:none)')
           elseif (maxval(fault%array(:,3)).lt.1000.0d0) then
  -            ! Depth units are meters
  +            ! All depths are small, using km?
               write(stderr,*) 'read_inputs: all fault depths are less than +1000 meters (using ', &
                               'gf_model ',trim(gf_model),')'
           endif
   
           ! Fault dimensions are in meters
           if (maxval(fault%array(:,6)).lt.100.0d0.and.maxval(fault%array(:,7)).lt.100.0d0) then
  +            ! All dimensions are small, using km?
               write(stderr,*) 'read_inputs: all fault dimensions are less than 100 meters ',&
                               ' using gf_model ',trim(gf_model)
           endif
   
  +    ! Point source faults in an elastic half-space
       elseif (gf_model.eq.'okada_pt') then
  +
  +        ! Depth is positive down, units in meters
           if (minval(fault%array(:,3)).lt.0.0d0) then
  -            ! Depth is positive down
  +            ! Found negative depth
               write(stderr,*) 'read_inputs: found fault depth of ',minval(fault%array(:,3))
  -            call usage('Depth is defined positive down for gf_model '//trim(gf_model))
  +            call usage('Depth is defined positive down for gf_model '//trim(gf_model)// &
  +                       '(usage:none)')
           elseif (maxval(fault%array(:,3)).lt.1000.0d0) then
  -            ! Depth units are meters
  +            ! All depths are small, using km?
               write(stderr,*) 'read_inputs: all fault depths are less than +1000 meters (using ', &
                               'gf_model ',trim(gf_model),')'
           endif
   
           ! Fault areas are in square meters
           if (maxval(fault%array(:,6)).lt.10000.0d0) then
  +            ! Area is small, using km^2?
               write(stderr,*) 'read_inputs: all fault areas are less than 100x100 square meters ',&
                               'using gf_model ',trim(gf_model)
           endif
   
  +    ! Triangular faults in an elastic half-space
       elseif (gf_model.eq.'triangle') then
           if (minval([fault%array(:,3),fault%array(:,6),fault%array(:,9)]).lt.0.0d0) then
  -            ! Depth is positive down
  +            ! Found negative depth
               write(stderr,*) 'read_inputs: found triangle vertex depth of ', &
                               minval([fault%array(:,3),fault%array(:,6),fault%array(:,9)])
  -            call usage('Depth is defined positive down for gf_model '//trim(gf_model))
  +            call usage('Depth is defined positive down for gf_model '//trim(gf_model)// &
  +                       '(usage:none)')
           endif
   
           if (maxval([fault%array(:,3),fault%array(:,6),fault%array(:,9)]).lt.1000.0d0) then
  -            ! Depth units are meters
  +            ! All depths are small, using km?
               write(stderr,*) 'read_inputs: all fault depths are less than +1000 meters (using ', &
                               'gf_model ',trim(gf_model),')'
           endif
  
positive.code:     endif
  else
      if (inversion_mode.eq.'anneal') then
          call usage('read_inputs: a slip constraint file is required for annealing search')
      elseif (inversion_mode.eq.'anneal-psc') then
          call usage('read_inputs: a slip constraint file is required for annealing+'//&
                     'pseudo-coupling search')
      endif
  endif
  
  
negative.code:         endif
          allocate(rake_constraint%array(fault%nrows,rake_constraint%ncols),stat=ierr)
          if (ierr.ne.0) then
              call usage('read_inputs: error allocating memory to rake_constraint%array (usage:none)')
          endif
          rake_constraint%nrows = fault%nrows
          rake_constraint%array(:,1) = dp1
  