ReadNetworkGroundwater Subroutine

SUBROUTINE ReadNetworkGroundwater &
!
( filename, domain, network )

!arguments with intent in:
CHARACTER (LEN = *), INTENT(IN) :: filename
TYPE(grid_integer), INTENT(IN)  :: domain

!arguments with intent out:
TYPE(NetworkGroundwater), INTENT(OUT) :: network


!local declarations
INTEGER  :: k
INTEGER  :: fileunit
INTEGER  :: err_io
INTEGER  :: id
LOGICAL  :: check_domain

!-------------------------------end of declaration----------------------------- 
!open file in readonly mode
fileunit = GetUnit ()
OPEN (unit = fileunit, file = filename, status = "old", &
     action = "read", iostat = err_io )

IF (err_io /= 0) THEN
  CALL Catch ('error', 'HydroNetwork',        &
              'error opening file of groundwater hydro network: ',    &
              code = openFileError, argument = filename )
END IF

!count number of branches in  network
k = 0
READ(fileunit,*)
READ(fileunit,*)

DO !loop till end of file	
  READ(fileunit,*,iostat=err_io) id
		IF (err_io /= 0) THEN
			EXIT
    ELSE
      k = k + 1
    END IF
END DO

IF (k < 1) THEN
    CALL Catch ('error', 'HydroNetwork', &
               'no branches detected in groundwater network') 
ELSE
   CALL Catch ('info', 'HydroNetwork', &
               ' number of detected branches in groundwater network: ', &
               argument = ToString(k) ) 
END IF

!rewind file
REWIND (fileunit)
!skip first two lines
READ(fileunit,*)
READ(fileunit,*)

!allocate branches
network % nreach = k
ALLOCATE (network % branch(k))

!read info from file
DO k = 1, network % nreach
  READ (fileunit,*) network % branch (k) % id, network % branch (k) % x0, &
                    network % branch (k) % y0, network % branch (k) % x1, &
                    network % branch (k) % y1, network % branch (k) % ncells, &
                    network % branch (k) % order, network % branch (k) % riverbed, &
                     network % branch (k) % scalefactor_conductivity 
  !compute i0, j0, i1, j1
  CALL GetIJ (X = network % branch (k) % x0, Y = network % branch (k) % y0, &
              grid = domain, i = network % branch (k) % i0, j = network % branch (k) % j0, &
              check = check_domain)
  IF (.NOT. check_domain) THEN
      CALL Catch ('error', 'HydroNetwork', &
                   'beginning cell out of domain in branch of groundwater with id: ', &
                   argument = ToString(network % branch(k) % id ) )
  END IF

  CALL GetIJ (X = network % branch (k) % x1, Y = network % branch (k) % y1, &
              grid = domain, i = network % branch (k) % i1, j = network % branch (k) % j1, &
              check = check_domain )
  IF (.NOT. check_domain) THEN
      CALL Catch ('error', 'HydroNetwork', &
                   'end cell out of domain in branch of groundwater with id: ', &
                   argument = ToString(network % branch(k) % id ) )
  END IF
END DO

!close file
CLOSE (fileunit)

RETURN
END SUBROUTINE ReadNetworkGroundwater