module pdsi_types
  implicit none
  private
  integer, parameter :: dp = selected_real_kind(15, 307)
  real(dp), parameter :: MISSING = -99.0_dp

  ! Define maximum number of periods (e.g., 52 weeks or 12 months)
  integer, parameter :: MAX_PERIODS = 52
  integer, parameter :: MAX_YEARS = 1000

  ! Derived type to hold linked list node equivalent data (in arrays)
  type :: LinkedList_t
    real(dp), allocatable :: data(:)
    integer :: size = 0
    integer :: capacity = 0
  contains
    procedure :: init => ll_init
    procedure :: insert => ll_insert
    procedure :: remove_head => ll_remove_head
    procedure :: is_empty => ll_is_empty
    procedure :: get_size => ll_get_size
    procedure :: to_array => ll_to_array
    procedure :: clear => ll_clear
  end type LinkedList_t

contains

  subroutine ll_init(this, initial_capacity)
    class(LinkedList_t), intent(inout) :: this
    integer, intent(in), optional :: initial_capacity
    integer :: n
    
    if (present(initial_capacity)) then
      n = initial_capacity
    else
      n = 1024
    end if
    allocate(this%data(n))
    this%capacity = n
    this%size = 0
  end subroutine ll_init

  subroutine ll_insert(this, x)
    class(LinkedList_t), intent(inout) :: this
    real(dp), intent(in) :: x
    if (this%size == this%capacity) then
      call extend_capacity(this)
    end if
    this%size = this%size + 1
    this%data(this%size) = x
  end subroutine ll_insert

  subroutine extend_capacity(this)
    class(LinkedList_t), intent(inout) :: this
    real(dp), allocatable :: tmp(:)
    integer :: new_capacity

    new_capacity = this%capacity * 2
    allocate(tmp(new_capacity))
    tmp(1:this%size) = this%data(1:this%size)
    deallocate(this%data)
    this%data = tmp
    this%capacity = new_capacity
  end subroutine extend_capacity

  function ll_remove_head(this) result(x)
    class(LinkedList_t), intent(inout) :: this
    real(dp) :: x
    integer :: i

    if (this%size == 0) then
      x = MISSING
      return
    end if

    x = this%data(1)
    do i = 2, this%size
       this%data(i-1) = this%data(i)
    end do
    this%size = this%size -1
  end function ll_remove_head

  logical function ll_is_empty(this)
    class(LinkedList_t), intent(in) :: this
    ll_is_empty = (this%size == 0)
  end function ll_is_empty

  integer function ll_get_size(this)
    class(LinkedList_t), intent(in) :: this
    ll_get_size = this%size
  end function ll_get_size

  function ll_to_array(this) result(arr)
    class(LinkedList_t), intent(in) :: this
    real(dp), allocatable :: arr(:)
    if (this%size > 0) then
      allocate(arr(this%size))
      arr = this%data(1:this%size)
    else
      allocate(arr(0))
    end if
  end function ll_to_array

  subroutine ll_clear(this)
    class(LinkedList_t), intent(inout) :: this
    if (allocated(this%data)) deallocate(this%data)
    this%capacity = 0
    this%size = 0
  end subroutine ll_clear

end module pdsi_types
module pdsi_io
  use pdsi_types
  implicit none
  private
  public :: read_temperature_file, read_precip_file, read_normal_temp, read_parameter_file

contains

  ! Reads temperature file (monthly or weekly) into a 2D array: years x periods
  ! Assumes file format: year integer, followed by n_periods real temperature values
  function read_temperature_file(filename, n_years, n_periods) result(temp_array)
    character(len=*), intent(in) :: filename
    integer, intent(in) :: n_years, n_periods
    real(dp), allocatable :: temp_array(:,:)

    integer :: i, j, year_read
    character(len=256) :: line
    integer :: ios
    integer :: unit=20

    allocate(temp_array(n_years,n_periods))
    temp_array = MISSING

    open(unit=unit, file=filename, status='old', action='read', iostat=ios)
    if (ios /= 0) then
      print *, 'Error: Cannot open temperature file: ', filename
      stop
    end if

    do i = 1, n_years
      read(unit, *, iostat=ios) year_read
      if (ios /= 0) exit
      do j = 1, n_periods
        read(unit, *, iostat=ios) temp_array(i, j)
        if (ios /= 0) exit
      end do
    end do

    close(unit)

  end function read_temperature_file


  ! Reads precipitation file with the same format as temperature
  function read_precip_file(filename, n_years, n_periods) result(precip_array)
    character(len=*), intent(in) :: filename
    integer, intent(in) :: n_years, n_periods
    real(dp), allocatable :: precip_array(:,:)

    integer :: i, j, year_read
    integer :: ios
    integer :: unit=21

    allocate(precip_array(n_years,n_periods))
    precip_array = MISSING

    open(unit=unit, file=filename, status='old', action='read', iostat=ios)
    if (ios /= 0) then
      print *, 'Error: Cannot open precipitation file: ', filename
      stop
    end if

    do i = 1, n_years
      read(unit, *, iostat=ios) year_read
      if (ios /= 0) exit
      do j = 1, n_periods
        read(unit, *, iostat=ios) precip_array(i, j)
        if (ios /= 0) exit
      end do
    end do

    close(unit)

  end function read_precip_file


  ! Reads normal temperature array (average monthly or weekly temps) from single line or multiple lines
  function read_normal_temp(filename, n_periods) result(norm_temp)
    character(len=*), intent(in) :: filename
    integer, intent(in) :: n_periods
    real(dp), allocatable :: norm_temp(:)

    integer :: i, ios, unit=22
    real(dp) :: temp_val

    allocate(norm_temp(n_periods))
    norm_temp = MISSING

    open(unit=unit, file=filename, status='old', action='read', iostat=ios)
    if (ios /= 0) then
      print *, 'Error: Cannot open normal temperature file: ', filename
      stop
    end if

    do i = 1, n_periods
      read(unit, *, iostat=ios) temp_val
      if (ios /= 0) exit
      norm_temp(i) = temp_val
    end do

    close(unit)

  end function read_normal_temp


  ! Reads parameter file containing AWC and TLA on one or two lines
  subroutine read_parameter_file(filename, AWC, TLA)
    character(len=*), intent(in) :: filename
    real(dp), intent(out) :: AWC, TLA

    integer :: ios, unit=23
    real(dp) :: val1, val2

    AWC = MISSING
    TLA = MISSING

    open(unit=unit, file=filename, status='old', action='read', iostat=ios)
    if (ios /= 0) then
      print *, 'Warning: Cannot open parameter file: ', filename
      stop
    end if

    ! Try to read two numbers (either in one line or two lines)
    read(unit, *, iostat=ios) val1, val2
    if (ios /= 0) then
      ! try reading one number per line if two numbers on one line fail
      rewind(unit)
      read(unit, *, iostat=ios) val1
      if (ios /= 0) then
        print *, 'Error reading parameter file: ', filename
        stop
      end if
      read(unit, *, iostat=ios) val2
      if (ios /= 0) then
        print *, 'Error reading parameter file: ', filename
        stop
      end if
    end if

    ! Assign values assuming first is AWC, second TLA
    AWC = val1
    TLA = val2

    close(unit)
  end subroutine read_parameter_file

end module pdsi_io