      SUBROUTINE werrinit(fcerr,im,jnp,mlev,uerr,verr)
      include 'processors.h'
c      REAL uerr (im, jnp, mlev), verr (im, jnp, mlev), fcerr (
c     &im, jnp, mlev)
      REAL uerr (:, :, :), verr (:, :, :), fcerr (
     &:, :, :)

      COMMON / cortabs / hcorfe (1000, 3, 20, 2)


!hpf$ distribute fcerr(*,block,*) onto PROC_ARRAY_1D
!hpf$ distribute uerr(*,block,*) onto PROC_ARRAY_1D
!hpf$ distribute verr(*,block,*) onto PROC_ARRAY_1D

      real,dimension(im,jnp):: sinwt,phiwt,latitude,f,geostwt,tropwt
!hpf$ distribute phiwt(*,block) onto PROC_ARRAY_1D
!hpf$ distribute sinwt(*,block) onto PROC_ARRAY_1D
!hpf$ distribute latitude(*,block) onto PROC_ARRAY_1D
!hpf$ distribute f(*,block) onto PROC_ARRAY_1D
!hpf$ distribute geostwt(*,block) onto PROC_ARRAY_1D
!hpf$ distribute tropwt(*,block) onto PROC_ARRAY_1D

      DATA GRAVITY / 9.81 /

      print *, 'before all'
      pi = 4. * atan(1.)
      drc = pi / 180.
      CALL getconh('GRADEI  ',gradei)
      werror = -hcorfe(1,3,1,2)
      sqrtwe = sqrt(werror)
      print *, 'before forall'
      forall(j=1:jnp)
     &      latitude(:,j) = (180. / float(jnp - 1)) * (j - 1) - 90.
c      do j=1,jnp
c           latitude(:,j) = (180. / float(jnp - 1)) * (j - 1) - 90.
c       end do
c      print *, latitude
      phiwt = abs(latitude) / 30.
      latitude = latitude * drc
      f = 2. * (7.29E-5) * sin(latitude)
      sinwt = phiwt * 90. * drc
      geostwt = phiwt * sinwt
      tropwt = 1. - phiwt
      WHERE((phiwt * 30.) .GT. 30.)
              geostwt = 1.
              tropwt = 1.
      ENDWHERE
c   Here is my fix to the bug in the sequential code, div by 0.
      where(f.eq.0)
         f=f+0.000001
      end where
      print *, f
      print *, 'before DO'
      DO 7004 k = 1 , mlev
            uerr(:,:,k) = abs(gradei / f) * sqrtwe * fcerr(:,:,k)
            uerr(:,:,k) = geostwt * uerr(:,:,k) + tropwt * 5.
            verr(:,:,k) = uerr(:,:,k)
 7004 CONTINUE
      PRINT *, ' Geostrophic wind errors computed. '
      RETURN
      END
