C
C
CGH 20 Aug 1999 added Wrho_t 
C
      subroutine statistics(iter)
         integer iter
         
	 include 'CASCO.DIM'
	 include 'errors.h'
	 include 'covariance.h'
	 include 'statistics.h' 
	 include 'wrho.h'
	 include 'deltabcs.h'

	 integer nobs, nr,n_eobs
	 real    rms1, rms0
	 real    rwms1, rwms0
	 real    rmse1, rmse0
	 real    rwmse1, rwmse0
c	 real    vel_cost
c	 real    bc_cost

	 integer nobservations, nelev_bcs,nelevobservations
	 real    rms_2vec, rwms_2vec,rms_1vec
C        real    rwms_1vec
	 real    dWDELTAd, xWRHOx

	 nobs = nobservations()
	 n_eobs = nelevobservations()
	 nr   = nelev_bcs()

	 write(*,*)'statistics: iter: nobs ', iter, nobs

         rms1 = rms_2vec(uerr_best,verr_best,nobs)
         rms0 = rms_2vec(uerror,verror,nobs)
         write(*,*)'ITER, rms1, rms0 ', iter, rms1,rms0
         
	 rmse1 = rms_1vec(eerr_best,n_eobs)
         rmse0 = rms_1vec(eerror,n_eobs)
         write(*,*)'ITER, rmse1, rmse0 ', iter, rmse1,rmse0

         rwms1 = rwms_2vec(uerr_best,verr_best,W_delta,nobs)
         rwms0 = rwms_2vec(uerror,verror,W_delta,nobs)
         write(*,*)'ITER, rwms1, rwms0 ', iter, rwms1,rwms0

         rwmse1 = rwms_2vec(eerr_best,eerr_best,W_eps,n_eobs)
         rwmse0 = rwms_2vec(eerror,eerror,W_eps,n_eobs)
         write(*,*)'ITER, rwmse1, rwmse0 ', iter, rwmse1,rwmse0

C save cost function at previous iteration
         if(iter.gt.1) then 
	    costm1 = costfunc
	    bc_costm1 = bc_cost
	    vel_costm1 = vel_cost
	    elev_costm1 = elev_cost
	 endif 

C update cost function variables
c costfunc 
	 vel_cost = dWDELTAd(uerr_best,W_delta,nobs)+ 
     &              dWDELTAd(verr_best,W_delta,nobs)
     
	 elev_cost = dWDELTAd(eerr_best,W_eps,n_eobs) 
	 
	 bc_cost  = xWRHOx(dbc_best,Wrho,Wrho_t,IQrho,JQrho,nr,ntbc)

	 costfunc =  vel_cost + elev_cost + bc_cost  

C root weighted mean square velocity error
 	 rwms_vel  =  sqrt(dWDELTAd(uerr_best,W_delta,nobs)
     &                   + dWDELTAd(verr_best,W_delta,nobs))
         rms_vel   =  rms1
 	 rwms_elev =  rwmse1
 	 rms_elev  =  rmse1
     
         
     
         if(iter.eq.0) then 
         write(*,*)'costly mistakes',  iter,rwms_elev,rms_elev
            write(*,*)' COST:  iter, costfunc, rwms_vel,  rms_vel,',
     &	    '   rwms_elev,   rms_elev'
            write(*,101)iter,costfunc,rwms_vel,rms_vel,
     &rwms_elev,rms_elev
     
     	    cost0 = costfunc
	 else 
         write(*,*)'costly mistakes',  iter,rwms_elev,rms_elev
	    write(*,*)' COST-: ',iter-.5, rwms0**2,rwms0,rms0  
            write(*,101)iter,costfunc,rwms_vel,rms_vel,
     &rwms_elev,rms_elev
            write(*,*)'COST+:',costfunc,vel_cost,elev_cost,bc_cost
	 endif
	 write(*,*)'COST-bcs  COST-vel COST-elev: ', bc_cost,vel_cost,
     &      elev_cost  

101      format(' COST: ', I5, 5(1x,G10.4))

      return 
      end
