!! Diagnostics: U & V on earth coordinates

MODULE module_calc_uvmet

  CONTAINS
  SUBROUTINE calc_uvmet(SCR4, cname, cdesc, cunits)

  USE constants_module
  USE module_model_basics

  !Arguments
  real,allocatable,dimension(:,:,:,:)             :: SCR4
  character (len=128)                             :: cname, cdesc, cunits

  !Local
  integer                                         :: i, j
  real                                            :: cone
  real, dimension(west_east_dim,south_north_dim)  :: diff, alpha

  
  SCR4 = 0.0


  IF ( map_proj .ge. 3 ) THEN                         ! No need to rotate
    SCR4(:,:,:,1) = UUU
    SCR4(:,:,:,2) = VVV
  END IF


  cone = 1.                                          !  PS
  IF ( map_proj .eq. 1) THEN                         !  Lambert Conformal mapping
    IF (ABS(truelat1-truelat2) .GT. 0.1) THEN
       cone=(ALOG(COS(truelat1*RAD_PER_DEG))-            &
             ALOG(COS(truelat2*RAD_PER_DEG))) /          &
       (ALOG(TAN((90.-ABS(truelat1))*RAD_PER_DEG*0.5 ))- &
        ALOG(TAN((90.-ABS(truelat2))*RAD_PER_DEG*0.5 )) )
    ELSE
       cone = SIN(ABS(truelat1)*RAD_PER_DEG )
    ENDIF
  END IF


  diff = XLONG - stand_lon
  DO i = 1, west_east_dim
  DO j = 1, south_north_dim
    IF ( diff(i,j) .gt. 180. ) THEN
      diff(i,j) = diff(i,j) - 360.
    END IF
    IF ( diff(i,j) .lt. -180. ) THEN
      diff(i,j) = diff(i,j) + 360.
    END IF
  END DO
  END DO


  DO i = 1, west_east_dim
  DO j = 1, south_north_dim
     IF ( XLAT(i,j) .lt. 0. ) THEN
       alpha(i,j) = - diff(i,j) * cone * RAD_PER_DEG
     ELSE
       alpha(i,j) = diff(i,j) * cone * RAD_PER_DEG
     END IF
  END DO
  END DO


  DO k = 1,bottom_top_dim
    SCR4(:,:,k,1) = VVV(:,:,k)*sin(alpha) + UUU(:,:,k)*cos(alpha)
    SCR4(:,:,k,2) = VVV(:,:,k)*cos(alpha) - UUU(:,:,k)*sin(alpha)
  END DO


  cname    = "uvmet"
  cdesc    = "Rotated wind component"
  cunits   = "m s-1"
  
  END SUBROUTINE calc_uvmet

END MODULE module_calc_uvmet