module gt3map

    use gt3read, only: gt3_file, gt3_header, open, close
    use dc_types, only: string, token
    use dc_trace, only: message
    use findfile, only: axis_t

implicit none

    type timetab_t
        double precision:: time
    end type


    ! 軸表のスキーマ
    type axistab_t
        ! 軸同定情報
        type(axis_t):: spec
        ! 軸変数の netCDF 名前
        character(len = token), pointer:: name
        ! 軸重み変数の netCDF 名前
        character(len = token), pointer:: name_w
        ! 軸位置、重みのヘッダ
        type(gt3_header):: header, header_w
        ! 軸位置、重みの値
        real, pointer:: value(:), value_w(:)
        ! 上下端の差 (周期的な場合だけ modulo として使用)
        real:: range
    end type

    type varspec_t
        ! 変数自体の要素名
        character(len = 16):: item
        ! 軸
        type(axis_t):: axis(3)
    end type

    type vartab_t
        character(len = token), pointer:: name
        type(gt3_header):: header
        type(varspec_t):: spec
    end type

    ! 元ファイルにあったすべてのユニットについてのデータベース
    type unitlist_t
        ! 記録番号
        integer:: recno
        ! ヘッダ部
        type(gt3_header):: header
        ! この記録を無視すべきか
        logical:: valid
        type(unitlist_t), pointer:: next
        type(unitlist_t), pointer:: prev
    end type

    ! ある GTOOL3 レコードストリームから gtool4 ファイルを構成するための情報
    ! これを使って直接読みだしも、理論的には可能
    type gt3_map
        type(gt3_file):: file
        ! 必須属性
        character(len = token):: title, source, institution
        character(len = string):: history
        ! 時刻の単位
        character(len = string):: time_unit
        ! ユニット・リスト
        type(unitlist_t), pointer:: unit
        integer:: n_units
        ! 変数リスト
        type(vartab_t), pointer:: vartab(:)
        integer:: n_vars
        ! 時刻リスト
        type(timetab_t), pointer:: timetab(:)
        integer:: n_times
        ! 軸リスト
        type(axistab_t), pointer:: axistab(:)
        integer:: n_axes
        ! 変数、時刻など gtool/NetCDF 変数の名前を一括管理する表
        character(len = token), pointer:: nametab(:)
    end type

    interface open
        module procedure gt3mapOpen
    end interface

    interface close
        module procedure gt3mapClose
    end interface

    interface get_axis
        module procedure gt3mapgetaxis
    end interface

    interface get_variable
        module procedure gt3mapgetvarinfo
    end interface

    interface getHeader
        module procedure gt3mapgetheader
    end interface

    interface operator(==)
        module procedure axisspec_eqv
        module procedure varspec_eqv
    end interface

    integer, parameter:: GT3_DUPITEM = 1000

contains

    logical function axisspec_eqv(a, b) result(result)
        type(axis_t), intent(in):: a, b
    continue
        result = (a%item == b%item) .and. (a%start == b%start) .and. &
            & (a%end == b%end)
    end function

    logical function varspec_eqv(a, b) result(result)
        type(varspec_t), intent(in):: a, b
    continue
        result = a%item == b%item
        result = result .and. all(a%axis(:)%item == b%axis(:)%item)
        result = result .and. all(a%axis(:)%start == b%axis(:)%start)
        result = result .and. all(a%axis(:)%end == b%axis(:)%end)
    end function

    function varspec_from_header(header) result(result)
        use gt3read, only: gt3_header
        type(gt3_header), intent(in):: header
        type(varspec_t):: result
        result%item = header%item
        result%axis(:)%item = header%axis_item(:)
        result%axis(:)%start = header%axis_start(:)
        result%axis(:)%end = header%axis_end(:)
    end function

    ! map 構造体に1個の GTOOL3 ヘッダを登録する。これは、
    ! 一時的に使用される両方向リンクトリスト map%unit に格納される
    subroutine StoreHeader(map, header, irec)
        use gt3read, only: gt3_header
        type(gt3_map), intent(inout):: map
        type(gt3_header), intent(in):: header
        integer, intent(in):: irec
        type(unitlist_t), pointer:: unit
    continue
        allocate(unit)
        if (.not. associated(map%unit)) then
            unit%next => unit
            unit%prev => unit
            map%unit => unit
        else if (associated(map%unit%next, map%unit)) then
            unit%next => map%unit
            unit%prev => map%unit
            map%unit%prev => unit
            map%unit%next => unit
        else
            unit%next => map%unit
            unit%prev => map%unit%prev
            map%unit%prev%next => unit
            map%unit%prev => unit
        endif
        map%n_units = map%n_units + 1
        unit%recno = irec
        unit%header = header
        unit%valid = .true.
    end subroutine

    ! 一時ユニット表 map%unit から時刻の表 map%timetab を作成する。
    ! 時刻の数は予期できないが、いったん両方向リンクトリストを
    ! 作成し、この大きさで結果配列 map%timetab を確保した後で
    ! リンクトリストを開放しながら配列に転写していく。
    !
    subroutine build_time_table(map, ios)
        type(gt3_map), intent(inout):: map
        integer, intent(out):: ios
        type(unitlist_t), pointer:: unit
        type timelist_t
            real:: time
            type(timelist_t), pointer:: next
            type(timelist_t), pointer:: prev
        end type
        type(timelist_t), pointer:: timelist
        integer:: i
    continue
        ios = 0
        map%n_times = 0
        nullify(timelist)
        unit => map%unit
        if (.not. associated(unit)) return
        do
            if (unit%valid) then
                call insert_time(timelist, unit%header%time, ios)
            endif
            unit => unit%next
            if (associated(unit, map%unit)) exit
        enddo

        if (associated(timelist)) then
            allocate(map%timetab(map%n_times))
            nullify(timelist%prev%next)
            i = 1
            do
                map%timetab(i)%time = timelist%time
                timelist => timelist%next
                if (.not. associated(timelist)) exit
                deallocate(timelist%prev)
                i = i + 1
            enddo
        else
            nullify(map%timetab)
        endif

    contains

        subroutine insert_time(head, time, ios)
            type(timelist_t), pointer:: head
            integer, intent(in):: time
            integer, intent(out):: ios
            type(timelist_t), pointer:: timeent, cursor
        continue
            ios = 0
            cursor => head
            if (associated(cursor)) then
                do
                    if (cursor%time == time) return
                    if (cursor%time > time) exit
                    cursor => cursor%next
                    if (associated(cursor, head)) exit
                enddo
            endif
            allocate(timeent)
            map%n_times = map%n_times + 1
            timeent%time = unit%header%time
            if (.not. associated(cursor)) then
                timeent%next => timeent
                timeent%prev => timeent
                head => timeent
            else if (associated(cursor%next, cursor)) then
                timeent%next => cursor
                timeent%prev => cursor
                cursor%prev => timeent
                cursor%next => timeent
            else
                timeent%next => cursor
                timeent%prev => cursor%prev
                cursor%prev%next => timeent
                cursor%prev => timeent
            endif
        end subroutine

    end subroutine build_time_table

    ! 一時ユニット表 map%unit から変数表を作成する。
    ! 手順は時間表を作成するのにほぼ同じ。
    subroutine build_var_table(map, ios)
        type(gt3_map), intent(inout):: map
        integer, intent(out):: ios
        type(unitlist_t), pointer:: unit
        type varlist_t
            type(varspec_t):: spec
            type(gt3_header):: header
            type(varlist_t), pointer:: next
            type(varlist_t), pointer:: prev
        end type
        type(varlist_t), pointer:: varlist
        integer:: i
    continue
        unit => map%unit
        ios = 0
        map%n_vars = 0
        if (.not. associated(unit)) return
        do
            call insert_variable(varlist, unit%header, ios)
            if (ios == GT3_DUPITEM) then
                unit%valid = .false.
                ios = 0
            endif
            unit => unit%next
            if (associated(unit, map%unit)) exit
        enddo

        if (associated(varlist)) then
            allocate(map%vartab(map%n_vars))
            nullify(varlist%prev%next)
            i = 1
            do
                map%vartab(i)%spec = varlist%spec
                map%vartab(i)%header = varlist%header
                if (.not. associated(varlist%next)) exit
                varlist => varlist%next
                deallocate(varlist%prev)
                i = i + 1
            enddo
            deallocate(varlist)
        else
            nullify(map%vartab)
        endif

        return
    contains

        subroutine insert_variable(head, header, ios)
            type(varlist_t), pointer:: head
            type(gt3_header), intent(in):: header
            type(varspec_t):: spec
            integer, intent(out):: ios
            type(varlist_t), pointer:: varent
        continue
            spec = varspec_from_header(header)
            varent => head
            if (associated(varent)) then
                do
                    if (varent%spec == spec) return
                    if (varent%spec%item == spec%item) then
                        ios = GT3_DUPITEM
                        return
                    endif
                    varent => varent%next
                    if (associated(varent, head)) exit
                enddo
            endif
            allocate(varent)
            map%n_vars = map%n_vars + 1
            varent%spec = spec
            varent%header = header
            if (.not. associated(head)) then
                varent%next => varent
                varent%prev => varent
                head => varent
            else if (associated(head%next, head)) then
                varent%next => head
                varent%prev => head
                head%prev => varent
                head%next => varent
            else
                varent%next => head
                varent%prev => head%prev
                head%prev%next => varent
                head%prev => varent
            endif
        end subroutine

    end subroutine build_var_table

    ! 変数表から軸の表を作成する。
    ! 最後の軸 map%axistab(map%n_axes) は時間軸である。
    subroutine build_axis_table(map, ios)
        use findfile, only: load_axis
        type(gt3_map), intent(inout):: map
        integer, intent(out):: ios
        type axislist_t
            type(axis_t):: axis
            type(axislist_t), pointer:: next
            type(axislist_t), pointer:: prev
        end type
        type(axislist_t), pointer:: axislist
        integer:: i, j
    continue
        call message("listing axes")
        ios = 0
        nullify(axislist)
        map%n_axes = 1
        do, j = 1, map%n_vars
            do, i = 1, 3
                if (map%vartab(j)%spec%axis(i)%item /= '') then
                    call insert_axis(axislist, map%vartab(j)%spec%axis(i), ios)
                endif
            enddo
        enddo

        call message("building axis table (%d entry)", i=(/map%n_axes/))
        allocate(map%axistab(map%n_axes))
        if (associated(axislist)) then
            nullify(axislist%prev%next)
            i = 1
            do
                map%axistab(i)%spec = axislist%axis
                if (.not. associated(axislist%next)) exit
                axislist => axislist%next
                deallocate(axislist%prev)
                i = i + 1
            enddo
            deallocate(axislist)
        endif
        map%axistab(map%n_axes)%spec = axis_t("time", 1, map%n_times)

        call message("collecting axis files")
        do, i = 1, map%n_axes - 1
            call load_axis(map%axistab(i)%spec, 'LOC', &
                & map%axistab(i)%header, map%axistab(i)%value, ios, &
                & fullrange=map%axistab(i)%range)
            if (ios /= 0) return
            call load_axis(map%axistab(i)%spec, 'WGT', &
                & map%axistab(i)%header_w, map%axistab(i)%value_w, ios)
            if (ios /= 0) then
                write(*, '(A)') "weight for axis " // &
                    & map%axistab(i)%spec%item // " fot found"
                ios = 0
            endif
        enddo
        return
    contains

        subroutine insert_axis(head, spec, ios)
            type(axislist_t), pointer:: head
            type(axis_t), intent(in):: spec
            integer, intent(out):: ios
            type(axislist_t), pointer:: axisent
        continue
            ios = 0
            axisent => head
            if (associated(axisent)) then
                do
                    if (axisent%axis == spec) return
                    axisent => axisent%next
                    if (associated(axisent, head)) exit
                enddo
            endif
            map%n_axes = map%n_axes + 1
            allocate(axisent)
            axisent%axis = spec
            if (.not. associated(head)) then
                axisent%next => axisent
                axisent%prev => axisent
                head => axisent
            else if (associated(head%next, head)) then
                axisent%next => head
                axisent%prev => head
                head%prev => axisent
                head%next => axisent
            else
                axisent%next => head
                axisent%prev => head%prev
                head%prev%next => axisent
                head%prev => axisent
            endif
        end subroutine insert_axis

    end subroutine build_axis_table

    subroutine build_name_table(map, ios)
        use regex, only: match
        use dc_string, only: tolower, toChar
    implicit none
        type(gt3_map), intent(inout):: map
        integer, intent(out):: ios
        integer:: i, length, istart, j, offset
        character(len = token):: try
        character(len = *), parameter:: namecharset &
            & = 'abcdefghijklmnopqrstuvwxyz_0123456789'
    continue
        ! 変数名表の確保
        ! (軸数を倍にしているのは平均重みが入るから)
        allocate(map%nametab(map%n_vars + map%n_axes * 2 - 1))
        ! データ変数名。初期値として ITEM を入れる
        do, i = 1, map%n_vars
            map%nametab(i) = map%vartab(i)%spec%item
            map%vartab(i)%name => map%nametab(i)
        enddo
        ! 軸変数名の初期化 (軸重み変数名はあとまわし)
        offset = map%n_vars
        do, i = 1, map%n_axes
            map%nametab(i + offset) = map%axistab(i)%spec%item
            map%axistab(i)%name => map%nametab(i + offset)
        enddo
        ! まず GTOOL3 慣例の特定軸名称を置換する。末尾数字は残す
        do, i = 1, map%n_vars + map%n_axes
            try = map%nametab(i)
            call match('^GLON#d+$', try, start=istart, length=length)
            if (istart > 0) then
                map%nametab(i) = 'lon' // try(5: )
            endif
            call match('^GGLA#d+$', try, start=istart, length=length)
            if (istart > 0) then
                map%nametab(i) = 'lat' // try(5: )
            endif
            call match('^CSIG#d+#.[PM]$', try, start=istart, length=length)
            if (istart > 0) then
                map%nametab(i) = 'sigma' // try(length:length) &
                    & // try(5: length-2)
            endif
            ! 大文字から小文字へ
            try = map%nametab(i)
            call tolower(try)
            ! 名前文字以外の '_' への置換
            do
                istart = verify(try, namecharset // ' ')
                if (istart == 0) exit
                try(istart:istart) = '_'
            enddo
            ! 全部下線になると具合がわるいので
            istart = verify(try, ' _')
            if (istart == 0) then
                try = 'axis'
            else if (istart > 1) then
                try = try(istart: )
            endif
            map%nametab(i) = try
        enddo
        ! ついで名前の一意性を確保しつつ、簡略化を図る
        LOOP_U: do, i = 1, map%n_vars + map%n_axes
            ! もし数字を除去して一意名前が得られれば、そうする
            try = map%nametab(i)
            call match('^#a+#d+$', try, start=istart, length=length)
            if (istart > 0) then
                call match('^#a+', try, start=istart, length=length)
                try = try(1: length)
                if (all(map%nametab(:) /= try)) then
                    map%nametab(i) = try
                    cycle LOOP_U
                endif
            endif
            ! 逆に一意名前が得られないばあい...
            try = map%nametab(i)
            if (count(map%nametab(:) == try) == 1) cycle LOOP_U
            ! 数字で終わっていないなら数字を付加
            ! 数字で終わっているなら _数字 を付加
            call match('#d$', map%nametab(i), start=istart, length=length)
            do, j = 1, huge(1)
                if (istart == 0) then
                    try = trim(map%nametab(i)) // toChar(j)
                else
                    try = trim(map%nametab(i)) // '_' // toChar(j)
                endif
                if (all(map%nametab(:) /= try)) then
                    map%nametab(i) = try
                    cycle LOOP_U
                endif
            enddo
            ! もはや策なし。
            ios = -1
            return
        enddo LOOP_U

        ! 軸重み変数名
        offset = map%n_vars + map%n_axes
        do, i = 1, map%n_axes - 1
            map%nametab(i + offset) = trim(map%nametab(i + map%n_vars)) // "_wgt"
            map%axistab(i)%name_w => map%nametab(i + offset)
            if (.not. associated(map%axistab(i)%value_w)) then
                nullify(map%axistab(i)%name_w)
            endif
        enddo
        ! 時間については軸重みを付加しない
        nullify(map%axistab(map%n_axes)%name_w)
        ios = 0
    end subroutine

    subroutine gt3mapOpen(map, filename, ios)
        use gt3read, only: gt3_header, END_OF_FILE, &
            & getHeader, Open, Close, Rewind, SkipRecord
        use dc_string, only: printf
    implicit none
        type(gt3_map), intent(out):: map
        character(len = *), intent(in):: filename
        integer, intent(out):: ios
        type(gt3_header):: header
        integer:: rec_no, ios2
    continue
        !
        ! 初期化。ファイルを開く。
        nullify(map%unit, map%vartab, map%timetab, map%axistab)
        nullify(map%nametab)
        map%n_units = 0
        map%n_times = 0
        map%n_vars = 0
        map%n_axes = 0
        call open(map%file, filename, ios)
        if (ios /= 0) return
        !
        ! ファイル開けたら、全記録を読んでヘッダ情報を map に格納。
        call message("prescanning...")
        rec_no = 0
        do
            call GetHeader(map%file, header, ios)
            if (ios == END_OF_FILE) exit
            if (ios /= 0) goto 900
            rec_no = rec_no + 1
            if (rec_no == 1) then
                map%title = header%dataset
                call printf(map%source, "GTOOL3 format %d", i=(/header%idfm/))
                map%institution = header%modify_user
                map%time_unit = header%time_unit
                map%history = history_attr(header)
            endif
            call StoreHeader(map, header, rec_no)
            ! ユニット・リストを作成する
            call SkipRecord(map%file, ios)
            if (ios /= 0) goto 900
        enddo
        call message("%d units found.", i=(/rec_no/))
        ! ファイルは巻き戻すだけで閉じない。開かないで済む。
        ! (ただし、複数ファイルをサポートするならそうはいかない)
        call rewind(map%file, ios)
        if (ios /= 0) goto 900
        
        ! まず変数の表を作成する
        call build_var_table(map, ios)
        if (ios /= 0) goto 900
        ! ついで時刻の表を作成する
        call build_time_table(map, ios)
        if (ios /= 0) goto 900
        ! ついで次元の表を作成する
        call build_axis_table(map, ios)
        if (ios /= 0) goto 900
        ! 変数・次元に名前をつける
        call build_name_table(map, ios)
        if (ios /= 0) goto 900

        return
    900 continue
        call close(map%file, ios2)
        return
    contains

        function history_attr(h) result(result)
            character(len = string):: result
            type(gt3_header), intent(in):: h
            integer:: i, ilast
        continue
            result = trim(timefmt(h%create_date)) // " " // &
                & trim(h%create_user) // "> GTOOL3 create" // achar(10)
            ilast = 0
            do, i = 1, 8
                if (h%edit(i) == "") exit
                ilast = i
            enddo
            do, i = 1, ilast - 1
                result = trim(result) // "unknown unknown> GTOOL3 " // &
                    & trim(h%edit(i)) // ' ' // trim(h%edit_title(i)) // &
                    & achar(10)
            enddo
            i = ilast
            if (ilast > 0) then
                result = trim(result) // trim(timefmt(h%modify_date)) // &
                    & ' ' // trim(h%modify_user) // '> GTOOL3 ' // &
                    & trim(h%edit(i)) // ' ' // trim(h%edit_title(i)) // &
                    & achar(10)
            endif
            result = trim(result) // "unknown unknown> gt3conv" // &
                & achar(10)
        end function history_attr

    end subroutine gt3mapOpen

    function timefmt(gt3time) result(result)
        character(len = 19):: result
        character(len = *), intent(in):: gt3time
        character(len = token):: buffer
        integer:: y, m, d, h, mi, s, ios, i
    continue
        buffer = adjustl(gt3time)
        read(buffer, '(i4,i2,i2,1x3(i2))', iostat=ios) &
            & y, m, d, h, mi, s
        if (ios /= 0) then
            result = "unknown"
        else
            write(result, "(i4,2('-',i2),'T',2(i2,':'),i2)") &
                & y, m, d, h, mi, s
            do, i = 1, len(result)
                if (result(i:i) == ' ') result(i:i) = '0'
                if (result(i:i) == '*') result(i:i) = 'X'
            enddo
        endif
    end function

    subroutine gt3mapClose(map, ios)
        type(gt3_map), intent(inout):: map
        integer, intent(out):: ios
    continue
        call close(map%file, ios)
    end subroutine

    ! 属性のためにヘッダ情報を取得
    subroutine gt3mapGetHeader(map, class, order, header, varname)
    implicit none
        type(gt3_map), intent(in):: map
        character(len = *), intent(in):: class
        integer, intent(in):: order
        type(gt3_header), intent(out):: header
        character(len = *), intent(out):: varname
    continue
        select case(class)
            case("AXI")
                header = map%axistab(order)%header
                varname = map%axistab(order)%name
            case("WGT")
                header = map%axistab(order)%header_w
                varname = map%axistab(order)%name_w
            case("VAR")
                header = map%vartab(order)%header
                varname = map%vartab(order)%name
        end select
    end subroutine

    ! gtool_history の座標軸情報を得る
    subroutine gt3mapGetAxis(map, axis)
        use gtool_history, only: gt_history_axis
        type(gt3_map), intent(in):: map
        type(gt_history_axis), intent(out):: axis(:)
        integer:: i, length
    continue
        do, i = 1, size(axis)
            if (i == map%n_axes) then
                length = 0
                axis(i) = gt_history_axis( &
                    & map%axistab(i)%name, &
                    & length, &
                    & "time", &
                    & map%time_unit, &
                    & "float")
            else
                length = map%axistab(i)%spec%end - &
                    & map%axistab(i)%spec%start + 1
                axis(i) = gt_history_axis( &
                    & map%axistab(i)%name, &
                    & length, &
                    & map%axistab(i)%header%title, &
                    & map%axistab(i)%header%unit, &
                    & "float")
            endif
        enddo
    end subroutine

    ! すぐに historyAddVariable で使える変数情報 varinfo を発信
    subroutine gt3mapGetVarinfo(map, class, order, varinfo)
        use gtool_history, only: gt_history_varinfo
    implicit none
        type(gt3_map), intent(in):: map
        character(len = *), intent(in):: class
        integer, intent(in):: order
        type(gt_history_varinfo), intent(out):: varinfo
        character(len = token), pointer:: dimsp(:)
        integer:: ndims, stat, i, j
    continue
        select case (class)
        case('WGT')
            if (.not. associated(map%axistab(order)%name_w)) goto 999
            allocate(dimsp(1), stat=stat)
            if (stat /= 0) goto 999
            dimsp(1) = map%axistab(order)%name
            varinfo = gt_history_varinfo( &
                & map%axistab(order)%name_w, &
                & dimsp, &
                & map%axistab(order)%header_w%title, &
                & map%axistab(order)%header_w%unit, &
                & "float" )
        case('VAR')
            ndims = count(map%vartab(order)%spec%axis(:)%item /= "")
            allocate(dimsp(ndims + 1), stat=stat)
            if (stat /= 0) goto 999
            i = 1
            do, j = 1, 3
                dimsp(i) = map%vartab(order)%spec%axis(j)%item
                if (dimsp(i) == "") cycle
                dimsp(i) = varname_3to4(map, dimsp(i))
                i = i + 1
            enddo
            dimsp(ndims + 1) = map%axistab(map%n_axes)%name
            varinfo = gt_history_varinfo( &
                & map%vartab(order)%name, &
                & dimsp, &
                & map%vartab(order)%header%title, &
                & map%vartab(order)%header%unit, &
                & "float" )
        end select
        return
    999 continue
        nullify(dimsp)
        varinfo = gt_history_varinfo("", dimsp, "", "", "")
        return
    end subroutine

    function varname_3to4(map, item) result(result)
        type(gt3_map), intent(in):: map
        character(len = token):: result
        character(len = *), intent(in):: item
        integer:: i
        do, i = 1, map%n_axes
            if (map%axistab(i)%spec%item == item) then
                result = map%axistab(i)%name
                return
            endif
        enddo
        do, i = 1, map%n_vars
            if (map%vartab(i)%spec%item == item) then
                result = map%vartab(I)%name
                return
            endif
        enddo
        result = ""
    end function

end module gt3map