#$=head1 NAME 
#$
#$invint2 - Inverse linear interpolation;
#$
#$=head1 SYNOPSIS
#$
#$C<call invint1(niter,coord,ord,o1,d1,mm,mmov,eps,aa,doprec)>
#$  
#$=head1 PARAMETERS        
#$  
#$=over 4 
#$  
#$=item niter - integer
#$  
#$      Number of itterations
#$
#$=item coord -  C<real(:)> 
#$     
#$      Coordinates
#$  
#$=item o1,d1 -  real
#$     
#$      First model position and sampling
#$
#$=item ord   -  C<real(:)>
#$
#$      Data values
#$
#$=item mm    -  C<real(:)>
#$
#$      Output model
#$
#$=item mmov  -  C<real(:,:)>
#$  
#$      Model movie
#$
#$=item eps   -  real
#$     
#$      Epsilon value if doing preconditioning
#$  
#$=item aa    -  type(filter)
#$
#$      Preconditioning operator
#$
#$=item prec  -  logical
#$  
#$      Whether or not to do preconditioning
#$  
#$
#$=back
#$
#$=head1 DESCRIPTION
#$
#$Perform inverse linear interpolation
#$
#$
#$=head1 SEE ALSO
#$
#$L<lint1>,L<helicon>,L<polydiv>,L<solver>
#$
#$=head1 LIBRARY
#$
#$B<geef90>
#$
#$=cut
#$

module invint2 {   # Inverse linear interpolation
  use lint1
  use helicon                           # regularized    by helix   filtering
  use polydiv                           # preconditioned by inverse filtering
  use cgstep_mod + solver_mod
contains
  subroutine invint1( niter, coord,ord, o1,d1, mm,mmov, eps, aa, doprec) {
    logical,                  intent( in)  :: doprec
    integer,                  intent( in)  :: niter
    real,                     intent( in)  :: o1, d1, eps
    real,    dimension( :),   intent( in)  :: ord  
    type( filter),            intent( in)  :: aa
    real,    dimension( :),   intent( out) :: mm             
    real,    dimension( :,:), intent( out) :: mmov         # model movie
    real,    dimension( :),   pointer      :: coord        # coordinate
    call lint1_init( o1, d1, coord)
    if( doprec) {                                          # preconditioning
       call polydiv_init( size(mm), aa)
       call solver_prec( lint1_lop, cgstep, niter = niter, x = mm, dat = ord,  
                    prec = polydiv_lop, nprec = size(mm), eps = eps, 
                    xmov = mmov)
       call polydiv_close()
    } else {                                               # regularization
       call helicon_init( aa)
       call solver_reg(  lint1_lop, cgstep, niter = niter, x = mm, dat = ord, 
                    reg = helicon_lop, nreg = size(mm), eps = eps,
                    xmov = mmov)
    }                                               
    call cgstep_close()
  }
}