! vi: set sw=4 ts=8: module gt3conv_dimtable use gtool use gt3read use gt3conv_vartable implicit none private public DisposeTables, Dimension public BuildDimensions, BuildTime, PutDimensionData, PutTimeData public DIMTABLE_ENTRY, TIMETABLE_ENTRY, dims_table, time_table type DIMTABLE_ENTRY type(NC_DIMENSION):: dim type(NC_VARIABLE):: var character(len = 30):: name character(len = 16):: item integer:: lower, upper real, pointer:: values(:) end type type(DIMTABLE_ENTRY), save, pointer:: dims_table(:) type TIMETABLE_ENTRY integer, pointer:: values(:) type(NC_DIMENSION):: dim type(NC_VARIABLE):: var end type type(TIMETABLE_ENTRY), save, pointer:: time_table(:) interface Dimension; module procedure LookupDimension; end interface contains type(NC_DIMENSION) function LookupDimension(item) result(result) character(len = 16), intent(in):: item integer:: idim continue do, idim = 1, size(dims_table) if (dims_table(idim)%item == item) then result = dims_table(idim)%dim; return endif enddo result = NC_DIMENSION_ERROR() end function subroutine DisposeTables integer:: i continue call DisposeVariableTable do, i = 1, size(dims_table) deallocate(dims_table(i)%values) enddo deallocate(dims_table) do, i = 1, size(time_table) deallocate(time_table(i)%values) enddo deallocate(time_table) nullify(time_table) end subroutine type(VSTRING) function AxisFilename(item, weight) result(result) character(len = 16), intent(in):: item logical, intent(in), optional:: weight character(len = *), parameter:: GTAXDIR_L = '/var/spool/gt3-dcl5/' character(len = *), parameter:: GTAXDIR_W = 'C:/DENNOU/GT3/' character(len = 8):: prefix logical:: exist continue prefix = 'GTAXLOC.' if (present(weight)) then if (weight) prefix = 'GTAXWGT.' endif result = prefix // trim(item) inquire(file=char(result), exist=exist) if (exist) return result = GTAXDIR_L // prefix // trim(item) inquire(file=char(result), exist=exist) if (exist) return result = GTAXDIR_W // prefix // trim(item) inquire(file=char(result), exist=exist) if (exist) return result = GtoolGetenv('GTAXDIR') // prefix // trim(item) inquire(file=char(result), exist=exist) end function subroutine PutDimensionData integer:: idim continue do, idim = 1, size(dims_table) call Put(dims_table(idim)%var, dims_table(idim)%values) enddo end subroutine subroutine ReadAxisFile(dim) type(DIMTABLE_ENTRY), intent(inout):: dim type(GT3_FILE):: axisfile type(GT3_HEADER):: header real, pointer:: buffer(:, :, :) integer:: lb, ub type(NC_ATTRIBUTE):: attr continue call Open(axisfile, char(AxisFilename(dim%item))) call Get(axisfile, header, buffer) call Close(axisfile) lb = lbound(buffer, 1) ub = ubound(buffer, 1) if (dim%lower < lb) then print *, 'gt3conv#ReadAxisFile: requested lower bound ', & & dim%lower, ' is less than axis file value ', lb stop endif if (dim%upper > ub) then print *, 'gt3conv#ReadAxisFile: requested upper bound ', & & dim%upper, ' is greater than axis file value ', ub stop endif allocate(dim%values(lb: ub)) dim%values(:) = buffer(lb: ub, lbound(buffer, 2), lbound(buffer, 3)) if (header%dataset(1: 1) == 'C') then attr = Attribute(dim%var, 'topology') attr = 'circular' call Dispose(attr) endif if (header%unit /= '') then attr = Attribute(dim%var, 'unit') ! fake attr = UnitString(header) call Dispose(attr) endif end subroutine type(VSTRING) function Dimname3to4(item, usedname) result(result) character(len = *), intent(in):: item character(len = *), intent(in):: usedname(:) character(len = len(item)):: buffer type(VSTRING):: name, grids, suffix integer:: a, b, c continue buffer = adjustl(item) a = scan(buffer, '0123456789') name = buffer(1: a-1) b = verify(buffer(a: ), '0123456789') if (b == 0) b = len(buffer) grids = buffer(a: b) suffix = trim(buffer(b + 1: )) if (name == 'CSIG' .and. suffix == '.M') then result = 'sigma_bound' else if (name == 'CSIG') then result = 'sigma' else if (name == 'GLON') then result = 'lon' else if (name == 'GGLA') then result = 'lat' else result = item endif if (all(usedname(:) /= char(result))) return result = result // grids if (all(usedname(:) /= char(result))) return name = result do, c = 1, huge(c) result = name // itos(c) if (all(usedname(:) /= char(result))) return enddo stop 'gt3conv#Dimname3to4' end function type(VSTRING) function UnitString(header) type(GT3_HEADER), intent(in):: header continue if (header%unit == 'deg') then if (header%item(1: 4) == 'GLON') then UnitString = 'degree_east' else UnitString = 'degree_north' endif else UnitString = header%unit endif end function subroutine BuildDimensions(file) type(NC_FILE), intent(inout):: file type(DIMTABLE_ENTRY), pointer:: dims_tmp(:) integer:: iv, ia, idim, num_dims continue allocate(dims_table(3)) num_dims = 0 ! すべての変数について軸 ITEM 名を検査し、重複しないものを登録 do, iv = 1, size(vars_table) do, ia = 1, 3 call StoreDimension(vars_table(iv)%header, ia) enddo enddo ! 次元表の大きさを最低限に縮小 if (size(dims_table) > num_dims) then dims_tmp => dims_table allocate(dims_table(1: num_dims)) dims_table(:) = dims_tmp(1: num_dims) deallocate(dims_tmp) endif ! 登録されたすべての次元を作成 do, idim = 1, num_dims call MakeDimension(file, dims_table(idim), & & dims_table(1: idim-1)%name) enddo contains subroutine MakeDimension(file, dim, usedname) type(NC_FILE), intent(inout):: file type(DIMTABLE_ENTRY), intent(inout):: dim character(len = *), intent(in):: usedname(:) integer:: length continue dim%name = Dimname3to4(dim%item, usedname) dim%upper = max(dim%upper, dim%lower) length = dim%upper - dim%lower + 1 dim%dim = Dimension(file, trim(dim%name), length) call NetcdfAssert dim%var = Variable(file, trim(dim%name), NF_FLOAT, (/dim%dim/)) call NetcdfAssert call ReadAxisFile(dim) end subroutine subroutine StoreDimension(header, ia) type(GT3_HEADER), intent(in):: header integer, intent(in):: ia character(len = 16):: item integer:: idim type(DIMTABLE_ENTRY), pointer:: dims_tmp(:) continue ! 空なら飛ばす item = vars_table(iv)%header%axis_item(ia) if (item == '') return ! 登録済みなら上下端だけ確認 do, idim = 1, num_dims if (item == dims_table(idim)%item) then dims_table(idim)%lower = min(dims_table(idim)%lower, & & header%axis_start(ia)) dims_table(idim)%upper = max(dims_table(idim)%upper, & & header%axis_end(ia)) endif enddo ! 未登録なら次元表を確保し書き込み num_dims = num_dims + 1 if (num_dims > size(dims_table)) then allocate(dims_tmp(size(dims_table) + 2)) dims_tmp(1: size(dims_table)) = dims_table(:) deallocate(dims_table) dims_table => dims_tmp endif dims_table(num_dims)%item = item dims_table(num_dims)%lower = header%axis_start(ia) dims_table(num_dims)%upper = header%axis_end(ia) end subroutine end subroutine BuildDimensions subroutine PutTimeData(file) type(NC_FILE), intent(in):: file integer:: it integer, pointer:: values(:) continue do, it = 1, size(time_table) values => time_table(it)%values if (.not. NetcdfPutInt(time_table(it)%var, & & values, size(values))) then call NetcdfAssert() endif enddo end subroutine subroutine BuildTime(file) type(NC_FILE), intent(inout):: file character(len = 16):: axisname integer:: iv, num_times continue allocate(time_table(1)) num_times = 0 do, iv = 1, size(vars_table) call StoreTime enddo time_table => time_table(1: num_times) call CreateTimeDimensions contains subroutine CreateTimeDimensions integer:: it continue ! 時間軸がないばあい if (num_times <= 0) then time_table(:)%dim = NC_DIMENSION_ERROR() return endif ! 時間軸が一つしかないばあい time_table(1)%dim = & & Dimension(file, 'time', size(time_table(1)%values)) time_table(1)%var = & & Variable(file, 'time', NF_INT, (/time_table(1)%dim/)) if (num_times == 1) return ! 時間軸が複数あるばあい do, it = 2, num_times axisname = NetcdfName('time', it) time_table(it)%dim = & & Dimension(file, axisname, size(time_table(1)%values)) time_table(it)%var = & & Variable(file, axisname, NF_INT, (/time_table(it)%dim/)) enddo end subroutine subroutine StoreTime integer:: it, len integer, pointer:: t(:) type(TIMETABLE_ENTRY), pointer:: tt(:) continue len = count(units_table(:)%var == iv) if (len < 1) stop 'gt3conv_v#StoreTime barf' allocate(t(len)) t = pack(units_table(:)%time, units_table(:)%var == iv) ! 既存時間軸の探索 do, it = 1, num_times if (int_array_equiv(t, time_table(it)%values)) then vars_table(iv)%time_ord = it deallocate(t) return endif enddo ! 新しい時間軸の登録 num_times = num_times + 1 len = size(time_table) if (num_times > len) then allocate(tt(len + 1)) tt(1: len) = time_table(1: len) deallocate(time_table) time_table => tt endif vars_table(iv)%time_ord = num_times time_table(num_times)%values => t end subroutine end subroutine logical function int_array_equiv(a, b) result(result) integer, intent(in):: a(:), b(:) continue if (size(a) /= size(b)) then result = .FALSE.; return endif result = all(a(:) == b(:)) end function end module