dir_spread.f

Go to the documentation of this file.

!                   *********************
                    SUBROUTINE dir_spread
!                   *********************
!
     &( dirspr, f, ndire , nf, npoin2)
!
!***********************************************************************
! TOMAWAC   V7P2                                   05/06/2016
!***********************************************************************
!
!brief    COMPUTES THE MEAN DIRECTIONAL SPREAD (=DIRECTIONAL
!+                WIDTH) S IN DEGREES.
!
!history  M. BENOIT
!+        28/12/95
!+        V1P1
!+   CREATED
!
!history  M. BENOIT
!+        05/07/96
!+        V1P2
!+   MODIFIED
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        13/07/2010
!+        V6P0
!+   Translation of French comments within the FORTRAN sources into
!+   English comments
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        21/08/2010
!+        V6P0
!+   Creation of DOXYGEN tags for automated documentation and
!+   cross-referencing of the FORTRAN sources
!
!history  G.MATTAROLO (EDF - LNHE)
!+        28/06/2011
!+        V6P1
!+   Translation of French names of the variables in argument
!
!history  S.E.BOURBAN (HRW)
!+        05/06/2016
!+        V7P2
!+   Name of the subroutine changed from SPREAD (an intrinsic Fortran
!+   routine) into DIR_SPREAD
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| COSMOY         |<--| WORK TABLE
!| DIRSPR         |<--| MEAN DIRECTIONAL SPREAD
!| F              |-->| VARIANCE DENSITY DIRECTIONAL SPECTRUM
!| NF             |-->| NUMBER OF FREQUENCIES
!| NDIRE          |-->| NUMBER OF DIRECTIONS
!| NPOIN2         |-->| NUMBER OF POINTS IN 2D MESH
!| SINMOY         |<--| WORK TABLE
!| TAUXC          |<--| WORK TABLE
!| TAUXE          |<--| WORK TABLE
!| TAUXS          |<--| WORK TABLE
!| VARIAN         |<--| WORK TABLE
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE declarations_tomawac, ONLY : deupi,gradeg, freq, dfreq, tailf,
     &                                 costet, sintet
!
      USE interface_tomawac, ex_dir_spread => dir_spread
      IMPLICIT NONE
!
!.....VARIABLES IN ARGUMENT
!     """"""""""""""""""""
      INTEGER, INTENT(IN)    :: NF    , NDIRE , NPOIN2
      DOUBLE PRECISION, INTENT(IN)    :: F(npoin2,ndire,nf)
      DOUBLE PRECISION, INTENT(INOUT) :: DIRSPR(npoin2)
!
!.....LOCAL VARIABLES
!     """""""""""""""""
      INTEGER  IP    , JP    , JF
      DOUBLE PRECISION AUXC  , AUXS  , DFDTET, DTETAR, AUXI
      DOUBLE PRECISION SEUIL , COEFT
      DOUBLE PRECISION SINMOY, COSMOY
      DOUBLE PRECISION VARIAN, TAUXE
      DOUBLE PRECISION TAUXS, TAUXC
!
!
      seuil=1.d-20
      dtetar=deupi/dble(ndire)
!
      DO ip=1,npoin2
        cosmoy=0.d0
        sinmoy=0.d0
        varian=0.d0
!-----C-------------------------------------------------------C
!-----C  SUMS UP THE DISCRETISED PART OF THE SPECTRUM         C
!-----C-------------------------------------------------------C
        DO jf=1,nf
!
          tauxc=0.d0
          tauxs=0.d0
          tauxe=0.d0
          dfdtet=dfreq(jf)*dtetar
          DO jp=1,ndire
            auxc=costet(jp)*dfdtet
            auxs=sintet(jp)*dfdtet
            tauxc=tauxc+f(ip,jp,jf)*auxc
            tauxs=tauxs+f(ip,jp,jf)*auxs
            tauxe=tauxe+f(ip,jp,jf)*dfdtet
          ENDDO
          cosmoy=cosmoy+tauxc
          sinmoy=sinmoy+tauxs
          varian=varian+tauxe
        ENDDO                    ! JF
!-----C-------------------------------------------------------------C
!-----C  TAKES INTO ACCOUNT THE HIGH FREQUENCY PART (OPTIONAL)      C
!-----C-------------------------------------------------------------C
        IF (tailf.GT.1.d0) THEN
          coeft=freq(nf)/((tailf-1.d0)*dfreq(nf))
          cosmoy=cosmoy+tauxc*coeft
          sinmoy=sinmoy+tauxs*coeft
          varian=varian+tauxe*coeft
        ENDIF
!-----C-------------------------------------------------------------C
!-----C  COMPUTES THE DIRECTIONAL WIDTH                             C
!-----C-------------------------------------------------------------C
        IF (varian.GT.seuil) THEN
          auxs=sinmoy/varian
          auxc=cosmoy/varian
          auxi=min(sqrt(auxs*auxs+auxc*auxc),1.d0)
          dirspr(ip)=sqrt(2.d0*(1.d0-auxi))*gradeg
        ELSE
          dirspr(ip)=seuil
        ENDIF
      ENDDO ! IP
!
      RETURN
      END