【メモ帳】gif 画像出力

gif 画像出力 module 分割

parameterized derived type は使わないことにした。よくよく考えると、元々 len 変数にサイズが依存する要素が無かったDeath！　

元々 parameterized derived type は、1987年以前の Fortran 8X で提案されていた機能で、オブジェクト指向とは別の発想で同様な機能拡張を実現しようとしたものだと思われます。その頃から２０年 Fortran 2003 でオブジェクト指向の要素が取り入れられた後では、配列長 len に関してはコンストラクタを定義しての allocatation が使え、種別 kind に関しては相応のクラス定義ないし万能クラス class(*) で対応できそうな気がします。

parameterized derived type は、プログラムも面倒な上に予想しにくい挙動もあり、gfortran ではまだ実装されていないようだし、Intel のコンパイラでもバグがなかなか取れないので、無理して使う必要もないかと思うこの頃。引数を持った型定義は、関数型言語などの影響か、流行っている気もしますが、Fortran では制限もあって使いにくいです。

出力結果

• 三角関数

• 八色タイル

Modern Fortran Explained (Numerical Mathematics and Scientific Computation)

ソース・プログラム

• 派生型の定義 　module m_gif_types

ユーザー定義型 I/O を使った。
function と interface でコンストラクタを定義して初期値を与えた。　

• LZW 圧縮エンコーダ　module m_gif_lzw

モジュール内での共有変数を使った。
アルゴリズムを簡素化。

• gif ファイル書き出し　module m_gif

将来アニメーション gif が出せるように image ブロック独立化の含みを少し持たせた。
アルゴリズムを簡素化。

    module m_gif_types
      use, intrinsic :: iso_fortran_env
      implicit none    
      interface t_gif_header
        procedure :: init_gif_header
      end interface 

      interface t_image_block
        procedure :: init_image_block
      end interface 

      type :: t_gif_header
        character(3)   :: signature = 'GIF'
        character(3)   :: version   = '89a' ! '87a'
        integer(int16) :: width  
        integer(int16) :: height 
        integer(int8)  :: pck = int(B'10001010', int8) ! 1:Global Col. 3:Color Res. 1:Sort Flg. 3:Size of Global Col.
        integer(int8)  :: background_color_index = 0
        integer(int8)  :: pixel_aspect_ratio     = 0
      contains
        procedure :: wr_gif_header
        generic   :: write(unformatted) => wr_gif_header
      end type  
       
      type :: t_image_block
        integer(int8)  :: image_separator = Z'2C'
        integer(int16) :: image_left_position = 00
        integer(int16) :: image_top_position  = 00
        integer(int16) :: image_width   
        integer(int16) :: image_height 
        integer(int8)  :: pck = int(B'00000000', int8)
        ! local color table
        !integer(int8) :: LZW_minimum_code_size 
        !integer(int8) :: block_size  
        !integer(int8) :: image_data
        !integer(int8) :: block_terminator = Z'00'
      contains
        procedure :: wr_image_block
        generic   :: write(unformatted) => wr_image_block
      end type t_image_block
      !
      private
      public :: t_gif_header, t_image_block
    contains
      type(t_gif_header) function init_gif_header(nx, ny, ngcol) result(res)
        integer, intent(in) :: nx, ny, ngcol
        res%width  = nx
        res%height = ny
        res%pck = iand(res%pck, B'11111000') + int(mod(ngcol, 8), int8) ! least 3bits : size of global color : 2**(ngcol + 1) 2..256
      end function init_gif_header
    
      type(t_image_block) function init_image_block(nx, ny) result(res)
        integer      , intent(in) :: nx, ny
        res%image_width  = nx
        res%image_height = ny
      end function init_image_block
    
      subroutine wr_gif_header(this, iw, iostat, iomsg)
        class(t_gif_header), intent(in) :: this
        integer, intent(in ) :: iw
        integer, intent(out) :: iostat
        character(*), intent(in out) :: iomsg
        write(iw) this%signature 
        write(iw) this%version   
        write(iw) this%width 
        write(iw) this%height
        write(iw) this%pck
        write(iw) this%background_color_index
        write(iw) this%pixel_aspect_ratio
        iostat = 0
      end subroutine wr_gif_header

      subroutine wr_image_block(this, iw, iostat, iomsg)
        class(t_image_block), intent(in) :: this
        integer, intent(in ) :: iw
        integer, intent(out) :: iostat
        character(*), intent(in out) :: iomsg
        write(iw) this%image_separator
        write(iw) this%image_left_position  
        write(iw) this%image_top_position 
        write(iw) this%image_width 
        write(iw) this%image_height 
        write(iw) this%pck 
        iostat = 0
      end subroutine wr_image_block
    end module m_gif_types 
    
    module m_gif_lzw
      use, intrinsic :: iso_fortran_env
      use m_gif_types
      implicit none
      type :: t_enc
        integer :: kbits = 0, id = 0
      end type    
      ! module variables
      integer :: dict(0:2**12 - 1), kd                ! 0-based,  used 0..kd
      integer :: nbits, kbits                         ! LZW minimum bit / LZW bit length
      !
      private
      public :: encoder
    contains
      subroutine clear_dict()
        integer :: i
        nbits = nbits
        kbits = nbits + 1
        kd    = 2**nbits + 1                          ! dictinary is used to kd
        dict  = 0
        forall(i = 0:kd) dict(i) = i                  ! dict(nbit) clear_code, dict(nbit+1) end_code 
      end subroutine clear_dict
    
      subroutine subenc(ienc, m1, m0)
        integer, intent(   out) :: ienc
        integer, intent(in    ) :: m1
        integer, intent(in out) :: m0
        integer :: k, id
        k = ishft(m0 + 1, 16) + m1                    ! ibuffer + 1 to avoid 0, 0 degeneracy 
        id = findloc(dict, k, dim = 1)                ! dictionary is 0-based 
        if (id == 0) then                             ! not found in the dictionary
          kd = kd + 1
          dict(kd) = k
          if (kd == 2**kbits + 1) kbits = kbits + 1
          ienc = m0
          m0   = m1
        else                                          ! found in the dictionary
          ienc = -1
          m0   = id - 1                               ! dictionary is 0-based
        end if    
      end subroutine subenc
    
      subroutine enc(irgb, code)
        integer, intent(in) :: irgb(:)
        type(t_enc), intent(out) :: code(:)
        integer :: i, ip, ic, ienc, ibuffer
        call clear_dict()                             ! clear dictionary and kd 
        ic = 1
        code(ic) = t_enc(kbits, 2**nbits)             ! clear_code 
        ibuffer = irgb(1)                             ! first pixel 
        do ip = 2, size(irgb)
          call subenc(ienc, irgb(ip), ibuffer) 
          if (ienc /= -1) then                        ! not found in the dictionary
            ic = ic + 1
            code(ic) = t_enc(kbits, ienc)
          end if  
          if (kd == 2**12 - 1) then                   ! dictionary full
            ic = ic + 1   
            code(ic) = t_enc(kbits, 2**nbits)         ! clear_code
            call clear_dict()                         ! clear dictionary 
          end if
        end do
        ic = ic + 1
        code(ic) = t_enc(kbits, ibuffer)              ! last pixel in buffer  
        ic = ic + 1 
        code(ic) = t_enc(kbits, 2**nbits + 1)         ! end_code  
      end subroutine enc
      
      subroutine encoder(nbs, irgb, icode)
        integer, intent(in ) :: nbs, irgb(:)
        integer(int8), allocatable, intent(out) :: icode(:)
        type(t_enc), allocatable :: code(:)
        integer :: i, j, k, ib, ic, nb
        nbits = nbs                                   ! nbits is a module variable   
        allocate(code(size(irgb) + 100))              ! ? large enough ?
        call enc(irgb, code)                          ! compress rgb to LZW code
        nb = ceiling(sum(code(:)%kbits) / 8.0)        !   required bytes  
        allocate(icode(nb))
        icode = 0                                     ! pack code to icode (bit string)    
        k = 0
        do i = 1, size(code)
          do j = 1, code(i)%kbits   
            if (btest(code(i)%id, j - 1)) then 
              ic = k / 8 + 1 
              ib = mod(k, 8)
              icode(ic) = ibset(icode(ic), ib)
            end if
            k = k + 1
          end do  
        end do
      end subroutine encoder 
    end module m_gif_lzw  
      
    module m_gif
      use, intrinsic :: iso_fortran_env
      use m_gif_types
      use m_gif_lzw
      implicit none
      private
      public :: wr_gif
    contains
      subroutine wr_gif(fn, irgb2)
        character(*), intent(in) :: fn
        integer     , intent(in) :: irgb2(:, :)
        type(t_gif_header ), allocatable :: head
        type(t_image_block), allocatable :: img    
        integer :: n_global_color 
        integer(int8), allocatable :: global_color_table(:, :)
        integer(int8), allocatable :: icode(:)
        integer :: nx, ny, nbits 
        ! set global color table
        n_global_color = 2
        allocate(global_color_table(3, 2**(n_global_color + 1)))
        global_color_table = 0
        global_color_table = reshape([[  0,  0,  0], [255,  0,  0], [0,255,  0], [  0,  0,255], &
                                      [255,255,  0], [255,  0,255], [0,255,255], [255,255,255]], [8,3])
        !  
        nx = size(irgb2, 1)
        ny = size(irgb2, 2)
        head = t_gif_header(nx, ny, n_global_color)  
        img  = t_image_block(nx, ny) 
        nbits = max(3, n_global_color + 1)           
        !
        block 
          integer :: ix, iy
          integer, allocatable :: irgb(:)
          allocate(irgb(size(irgb2)))
          forall (ix = 1:nx, iy = 1:ny) irgb(ix + nx * (iy - 1)) = irgb2(ix, iy)
          call encoder(nbits, irgb, icode)             ! rgb color to GIF LZW code
        end block  
        !
        block
          integer :: iw, i, nbyte
          open(newunit = iw, file = fn, access = 'stream')
          write(iw) head      
          write(iw) global_color_table
          write(iw) img     
          write(iw) achar(nbits)                     ! LZW_minimum_code_size = 2   (2+1=3bit)
          do i = 1, size(icode), 255         
            nbyte = min(i + 255, size(icode)) - i
            write(iw) achar(nbyte)             
            write(iw) icode(i:i+nbyte-1)
          end do
          write(iw) achar(00)                          ! block_terminator
          write(iw) achar(Z'3B')                     ! end of gif file
          close(iw)
        end block
      end subroutine wr_gif  
    end module m_gif
    
    program test
      use m_gif_types
      use m_gif
      implicit none
      integer, parameter :: nx = 512, ny = 512
      integer :: irgb(nx, ny)
      !
      block                       ! sin cos curve
        real :: x, y(5), pi = 4*atan(1.0)
        integer :: i
        do i = 1, 512
          x = 4 * pi / nx * i
          y = 255 * [sin(x), cos(x), sin(x / 2), cos(x / 2), 0.0] + 256
          irgb(i, int(y)) = [1, 2, 4, 6, 7]
        end do  
      end block
      call wr_gif('sincos.gif', irgb)
!      
      block                        ! tile  
        integer :: ix, iy
        forall (ix = 1:nx, iy = 1:ny) irgb(ix, iy) = mod(ix / 40 + iy / 40, 8) 
      end block
      call wr_gif('tile.gif', irgb)
    end program test

(H30.5.8 バグ修正)