New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 5581 for branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/OPA_SRC/OBS/obs_write.F90 – NEMO

Ignore:
Timestamp:
2015-07-10T13:28:53+02:00 (9 years ago)
Author:
timgraham
Message:

Merged head of trunk into branch

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/OPA_SRC/OBS/obs_write.F90

    r2287 r5581  
    66 
    77   !!---------------------------------------------------------------------- 
    8    !!   cdf_wri_p3d   : Write profile observation diagnostics in NetCDF format 
     8   !!   obs_wri_p3d   : Write profile observation diagnostics in NetCDF format 
    99   !!   obs_wri_sla   : Write SLA observation related diagnostics 
    1010   !!   obs_wri_sst   : Write SST observation related diagnostics 
    1111   !!   obs_wri_seaice: Write seaice observation related diagnostics 
    12    !!   cdf_wri_vel   : Write velocity observation diagnostics in NetCDF format 
     12   !!   obs_wri_vel   : Write velocity observation diagnostics in NetCDF format 
     13   !!   obs_wri_stats : Print basic statistics on the data being written out 
    1314   !!---------------------------------------------------------------------- 
    1415 
     
    3132   USE obs_sla_types 
    3233   USE obs_rot_vel          ! Rotation of velocities 
     34   USE obs_mpp              ! MPP support routines for observation diagnostics 
     35   USE lib_mpp        ! MPP routines 
    3336 
    3437   IMPLICIT NONE 
     
    256259      ! Write the obfbdata structure 
    257260      CALL write_obfbdata( cfname, fbdata ) 
    258        
     261 
     262      ! Output some basic statistics 
     263      CALL obs_wri_stats( fbdata ) 
     264 
    259265      CALL dealloc_obfbdata( fbdata ) 
    260266      
     
    414420      CALL write_obfbdata( cfname, fbdata ) 
    415421 
     422      ! Output some basic statistics 
     423      CALL obs_wri_stats( fbdata ) 
     424 
    416425      CALL dealloc_obfbdata( fbdata ) 
    417426 
     
    565574      CALL write_obfbdata( cfname, fbdata ) 
    566575 
     576      ! Output some basic statistics 
     577      CALL obs_wri_stats( fbdata ) 
     578 
    567579      CALL dealloc_obfbdata( fbdata ) 
    568580 
     
    715727      CALL write_obfbdata( cfname, fbdata ) 
    716728 
     729      ! Output some basic statistics 
     730      CALL obs_wri_stats( fbdata ) 
     731 
    717732      CALL dealloc_obfbdata( fbdata ) 
    718733 
     
    722737      !!----------------------------------------------------------------------- 
    723738      !! 
    724       !!                     *** ROUTINE obs_wri_p3d  *** 
     739      !!                     *** ROUTINE obs_wri_vel  *** 
    725740      !! 
    726741      !! ** Purpose : Write current (profile) observation  
     
    903918      CALL write_obfbdata( cfname, fbdata ) 
    904919       
     920      ! Output some basic statistics 
     921      CALL obs_wri_stats( fbdata ) 
     922 
    905923      CALL dealloc_obfbdata( fbdata ) 
    906924      
     
    912930   END SUBROUTINE obs_wri_vel 
    913931 
     932   SUBROUTINE obs_wri_stats( fbdata ) 
     933      !!----------------------------------------------------------------------- 
     934      !! 
     935      !!                     *** ROUTINE obs_wri_stats  *** 
     936      !! 
     937      !! ** Purpose : Output some basic statistics of the data being written out 
     938      !! 
     939      !! ** Method  : 
     940      !!  
     941      !! ** Action  : 
     942      !! 
     943      !!      ! 2014-08  (D. Lea) Initial version  
     944      !!----------------------------------------------------------------------- 
     945 
     946      !! * Arguments 
     947      TYPE(obfbdata) :: fbdata 
     948 
     949      !! * Local declarations 
     950      INTEGER :: jvar 
     951      INTEGER :: jo 
     952      INTEGER :: jk 
     953 
     954!      INTEGER :: nlev 
     955!      INTEGER :: nlevmpp 
     956!      INTEGER :: nobsmpp 
     957      INTEGER :: numgoodobs 
     958      INTEGER :: numgoodobsmpp 
     959      REAL(wp) :: zsumx 
     960      REAL(wp) :: zsumx2 
     961      REAL(wp) :: zomb 
     962 
     963      IF (lwp) THEN 
     964         WRITE(numout,*) '' 
     965         WRITE(numout,*) 'obs_wri_stats :' 
     966         WRITE(numout,*) '~~~~~~~~~~~~~~~'  
     967      ENDIF 
     968 
     969      DO jvar = 1, fbdata%nvar 
     970         zsumx=0.0_wp 
     971         zsumx2=0.0_wp 
     972         numgoodobs=0 
     973         DO jo = 1, fbdata%nobs 
     974            DO jk = 1, fbdata%nlev 
     975               IF ( ( fbdata%pob(jk,jo,jvar) < 9999.0 ) .AND. & 
     976                  & ( fbdata%pdep(jk,jo) < 9999.0 ) .AND. & 
     977                  & ( fbdata%padd(jk,jo,1,jvar) < 9999.0 ) ) THEN 
     978        
     979             zomb=fbdata%pob(jk, jo, jvar)-fbdata%padd(jk, jo, 1, jvar) 
     980                  zsumx=zsumx+zomb 
     981                  zsumx2=zsumx2+zomb**2 
     982                  numgoodobs=numgoodobs+1 
     983          ENDIF 
     984            ENDDO 
     985         ENDDO 
     986 
     987         CALL obs_mpp_sum_integer( numgoodobs, numgoodobsmpp ) 
     988         CALL mpp_sum(zsumx) 
     989         CALL mpp_sum(zsumx2) 
     990 
     991         IF (lwp) THEN 
     992       WRITE(numout,*) 'Type: ',fbdata%cname(jvar),'  Total number of good observations: ',numgoodobsmpp  
     993       WRITE(numout,*) 'Overall mean obs minus model of the good observations: ',zsumx/numgoodobsmpp 
     994            WRITE(numout,*) 'Overall RMS obs minus model of the good observations: ',sqrt( zsumx2/numgoodobsmpp ) 
     995       WRITE(numout,*) '' 
     996         ENDIF 
     997  
     998      ENDDO 
     999 
     1000   END SUBROUTINE obs_wri_stats 
     1001 
    9141002END MODULE obs_write 
Note: See TracChangeset for help on using the changeset viewer.