!!  Copyright (C) Stichting Deltares, 2005-2019.
!!
!!  This file is part of iMOD.
!!
!!  This program is free software: you can redistribute it and/or modify
!!  it under the terms of the GNU General Public License as published by
!!  the Free Software Foundation, either version 3 of the License, or
!!  (at your option) any later version.
!!
!!  This program is distributed in the hope that it will be useful,
!!  but WITHOUT ANY WARRANTY; without even the implied warranty of
!!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!!  GNU General Public License for more details.
!!
!!  You should have received a copy of the GNU General Public License
!!  along with this program.  If not, see <http://www.gnu.org/licenses/>.
!!
!!  Contact: imod.support@deltares.nl
!!  Stichting Deltares
!!  P.O. Box 177
!!  2600 MH Delft, The Netherlands.
!!
MODULE MOD_UDF_UTL

 USE WINTERACTER
 USE RESOURCE
 USE MOD_DBL
 USE MOD_IDF_PAR, ONLY : IDFOBJ
 USE KDTREE2_MODULE  
 USE MESHMOD
 USE MODPLOT, ONLY : MPW,LEGENDOBJ
 USE MOD_UTL, ONLY: UTL_GETUNIT,UTL_TIMING,UTL_READINITFILE,UTL_IDFGETCLASS
   
 IMPLICIT NONE

 TYPE(TMESH), POINTER :: UDFMESH
 TYPE(KDTREE2), POINTER :: UDFTREE  
 TYPE(KDTREE2_RESULT), DIMENSION(1) :: UDFRES

CONTAINS

 !##==============================================================
 SUBROUTINE UDF_DEALLOCATEMESH()
 !##==============================================================
 IMPLICIT NONE
 INTEGER :: I
 
 IF(ASSOCIATED(UDFMESH))THEN
  IF(ASSOCIATED(UDFMESH%TAG)) DEALLOCATE(UDFMESH%TAG)
  IF(ASSOCIATED(UDFMESH%NODE))DEALLOCATE(UDFMESH%NODE)
  IF(ASSOCIATED(UDFMESH%ELEM))THEN
   DO I=1,SIZE(UDFMESH%ELEM)
    IF(ASSOCIATED(UDFMESH%ELEM(I)%TAG)) DEALLOCATE(UDFMESH%ELEM(I)%TAG)
    IF(ASSOCIATED(UDFMESH%ELEM(I)%TAG)) DEALLOCATE(UDFMESH%ELEM(I)%TAG)
    IF(ASSOCIATED(UDFMESH%ELEM(I)%ITAG))DEALLOCATE(UDFMESH%ELEM(I)%ITAG)
   ENDDO
   DEALLOCATE(UDFMESH%ELEM)
  ENDIF
  
  IF(ASSOCIATED(UDFMESH%NODE2ELEM))THEN
   DO I=1,SIZE(UDFMESH%NODE2ELEM)
    IF(ASSOCIATED(UDFMESH%NODE2ELEM(I)%NUM)) DEALLOCATE(UDFMESH%NODE2ELEM(I)%NUM)
   ENDDO
   DEALLOCATE(UDFMESH%NODE2ELEM)
  ENDIF
 ENDIF
 
 END SUBROUTINE UDF_DEALLOCATEMESH
 
 !##==============================================================
 LOGICAL FUNCTION UDF_OPEN(IDF,UDFFILE,IREAD,IU)
 !##==============================================================
 IMPLICIT NONE
 TYPE(IDFOBJ),INTENT(INOUT) :: IDF
 CHARACTER(LEN=*),INTENT(IN) :: UDFFILE
 INTEGER,INTENT(IN) :: IREAD
 INTEGER,INTENT(OUT) :: IU
 INTEGER :: IOS
 CHARACTER(LEN=256) :: GMSHFILE,DATFILE,LINE,DIR
    
 UDF_OPEN=.FALSE.
 
 IU=UTL_GETUNIT(); OPEN(IU,FILE=UDFFILE,ACTION='READ',IOSTAT=IOS)
 IF(IOS.NE.0)THEN
  CALL WMESSAGEBOX(OKONLY,EXCLAMATIONICON,COMMONOK,'Cannot read file '//TRIM(UDFFILE),'Error'); RETURN
 ENDIF

 DIR=UDFFILE(:INDEX(UDFFILE,'\',.TRUE.))

 IF(.NOT.UTL_READINITFILE('MSHFILE',LINE,IU,0))RETURN; READ(LINE,*) GMSHFILE

 IF(.NOT.UTL_READINITFILE('XMIN',LINE,IU,0))RETURN; READ(LINE,*) IDF%XMIN
 IF(.NOT.UTL_READINITFILE('XMAX',LINE,IU,0))RETURN; READ(LINE,*) IDF%XMAX
 IF(.NOT.UTL_READINITFILE('YMIN',LINE,IU,0))RETURN; READ(LINE,*) IDF%YMIN
 IF(.NOT.UTL_READINITFILE('YMAX',LINE,IU,0))RETURN; READ(LINE,*) IDF%YMAX

 IF(.NOT.UTL_READINITFILE('DATFILE',LINE,IU,0))RETURN; READ(LINE,*) DATFILE
 IF(.NOT.UTL_READINITFILE('ZMIN',LINE,IU,0))RETURN; READ(LINE,*) IDF%DMIN
 IF(.NOT.UTL_READINITFILE('ZMAX',LINE,IU,0))RETURN; READ(LINE,*) IDF%DMAX

 CLOSE(IU); IU=0

 IF(IREAD.GT.0)THEN

  GMSHFILE=TRIM(DIR)//'\'//TRIM(GMSHFILE)
  DATFILE =TRIM(DIR)//'\'//TRIM(DATFILE )
  
  !## read the gmshfile
  IF(.NOT.READGMSH(UDFMESH,GMSHFILE,IREAD))RETURN

  UDFMESH%XMIN=IDF%XMIN
  UDFMESH%YMIN=IDF%YMIN
  UDFMESH%XMAX=IDF%XMAX
  UDFMESH%YMAX=IDF%YMAX

  !## read data and create kd-tree
  IF(IREAD.EQ.2)THEN

  !## create kd-tree
! CALL DATE_AND_TIME(VALUES=T1)    
   UDFTREE=>KDTREE2_CREATE(UDFMESH%NODE,SORT=.FALSE.,REARRANGE=.FALSE.) 
! CALL DATE_AND_TIME(VALUES=T2)
! CALL UTL_TIMING(T1,T2,TD)

  ENDIF

  !## open the datfile
  IU=UTL_GETUNIT(); OPEN(UNIT=IU,FILE=DATFILE,FORM='UNFORMATTED',STATUS='OLD',ACCESS='STREAM',IOSTAT=IOS)
  IF(IOS.NE.0)THEN
   CALL WMESSAGEBOX(OKONLY,EXCLAMATIONICON,COMMONOK,'Cannot read file '//TRIM(UDFFILE),'Error'); RETURN
  ENDIF
     
! !## get data
! IF(IREAD.EQ.0)THEN
!  IDF%DMIN= HUGE(1.0D0)
!  IDF%DMAX=-HUGE(1.0D0)
!  DO I=1,UDFMESH%NELEM
!   READ(IU) X
!   IDF%DMIN = MIN(IDF%DMIN,X)
!   IDF%DMAX = MAX(IDF%DMAX,X)
!  ENDDO
!  CLOSE(IU)
! ENDIF

 ENDIF
      
 UDF_OPEN=.TRUE.

 END FUNCTION UDF_OPEN
  
 !##==============================================================
 SUBROUTINE UDF_PLOTNETWORK(IU,LEG,NODATA,UMIN,UMAX)
 !##==============================================================
 IMPLICIT NONE
 INTEGER,INTENT(IN) :: IU
 REAL(KIND=DP_KIND),INTENT(IN) :: NODATA
 REAL(KIND=DP_KIND),INTENT(OUT) :: UMIN,UMAX
 TYPE(LEGENDOBJ) :: LEG
 INTEGER :: I,J,K,IRASTER
 REAL(KIND=DP_KIND) :: VAL
 REAL(KIND=DP_KIND),DIMENSION(:),ALLOCATABLE :: XP,YP
 INTEGER :: N,ICLR
 
 CALL WINDOWSELECT(0); IRASTER=WMENUGETSTATE(ID_IDFRASTERLINES,2)

 ALLOCATE(XP(4),YP(4))
 
! IF(IRASTER.EQ.0)CALL IGRFILLPATTERN(OUTLINE)
! IF(IRASTER.EQ.1)CALL IGRFILLPATTERN(SOLID)
 CALL IGRCOLOURN(WRGB(0,0,0))
 
 DO I=1,UDFMESH%NELEM
  N=UDFMESH%ELEM(I)%N
  DO J=1,N
   K    =UDFMESH%ELEM(I)%NUM(J)
   XP(J)=UDFMESH%NODE(1,K) 
   YP(J)=UDFMESH%NODE(2,K) 
  ENDDO
  IF(MAXVAL(XP(1:N)).GT.MPW%XMIN.AND.MINVAL(XP(1:N)).LT.MPW%XMAX.AND. &
     MAXVAL(YP(1:N)).GT.MPW%YMIN.AND.MINVAL(YP(1:N)).LT.MPW%YMAX)THEN  
   CALL IGRFILLPATTERN(SOLID)
   J=4*(I-1)+1
   READ(IU,POS=J) VAL
   IF(VAL.NE.NODATA.AND.VAL.EQ.VAL)THEN
    ICLR=UTL_IDFGETCLASS(LEG,VAL)
    UMIN=MIN(UMIN,VAL)
    UMAX=MAX(UMAX,VAL)
    CALL IGRCOLOURN(ICLR)
   ELSE
    !## white for nodata
    CALL IGRCOLOURN(WRGB(255,255,255))
   ENDIF
   CALL DBL_IGRPOLYGONCOMPLEX(XP,YP,N)
   IF(IRASTER.EQ.1)THEN
    CALL IGRCOLOURN(WRGB(0,0,0))
    CALL IGRFILLPATTERN(OUTLINE)
    CALL DBL_IGRPOLYGONCOMPLEX(XP,YP,N)
   ENDIF
  ENDIF

! IF (LFOUND) THEN
!  I = 4*(ELEMFOUND-1)+1 ! BYTES
!  READ(IU,POS=I) VAL
! ELSE    
!  VAL = NODATA
! END IF   

  
!  CALL DBL_IGRMOVETO(X1,Y1)
!  DO J=2,UDFMESH%ELEM(I)%N
!   K=UDFMESH%ELEM(I)%NUM(J)
!   X=UDFMESH%NODE(1,K) 
!   Y=UDFMESH%NODE(2,K) 
!   CALL DBL_IGRLINETO(X,Y)
!  ENDDO
!  CALL DBL_IGRLINETO(X1,Y1)
 ENDDO
 
 DEALLOCATE(XP,YP)
   
 END SUBROUTINE UDF_PLOTNETWORK
 
 !##==============================================================
 SUBROUTINE READUDFXY(X,Y,NODATA,VAL,IU)
 !##==============================================================
 IMPLICIT NONE
 REAL(KIND=DP_KIND), INTENT(IN) :: X
 REAL(KIND=DP_KIND), INTENT(IN) :: Y
 REAL(KIND=DP_KIND), INTENT(IN) :: NODATA
 REAL(KIND=DP_KIND), INTENT(OUT) :: VAL
 INTEGER,INTENT(IN) :: IU
 INTEGER :: NOD, ELEM, I, J, N, M, ELEMFOUND
 REAL(KDKIND), DIMENSION(2) :: QV
 LOGICAL :: LFOUND
 DOUBLEPRECISION, DIMENSION(2,4) :: V
 DOUBLEPRECISION, DIMENSION(2) :: P
  
 IF (X < UDFMESH%XMIN .OR. X > UDFMESH%XMAX .OR. Y < UDFMESH%YMIN .OR. Y > UDFMESH%YMAX) THEN
  VAL = NODATA
  RETURN  
 END IF
  
 QV(1) = X; QV(2) = Y
 P(1) = X; P(2) = Y
 !## find nearest node
 CALL KDTREE2_N_NEAREST(TP=UDFTREE,QV=QV,NN=1,RESULTS=UDFRES)
 NOD = UDFRES(1)%IDX
 !## loop over connected elements
 ELEMFOUND = 0
 LFOUND = .FALSE.
 DO I = 1, UDFMESH%NODE2ELEM(NOD)%N
  ELEM = UDFMESH%NODE2ELEM(NOD)%NUM(I)  
  N = UDFMESH%ELEM(ELEM)%N
  DO J = 1, N
   M = UDFMESH%ELEM(ELEM)%NUM(J)  
   V(1,J) = UDFMESH%NODE(1,M)
   V(2,J) = UDFMESH%NODE(2,M)
  END DO
  CALL POLYGON_CONTAINS_POINT_2D(N,V,P,LFOUND)
  IF (LFOUND) THEN
   ELEMFOUND = ELEM
   EXIT
  END IF  
 END DO    
  
 IF (LFOUND) THEN
  I = 4*(ELEMFOUND-1)+1 ! BYTES
  READ(IU,POS=I) VAL
 ELSE    
  VAL = NODATA
 END IF   
  
 END SUBROUTINE READUDFXY
  
END MODULE MOD_UDF_UTL