anchor.diff:        endif
       end select
     !---------------------------------------------------------------------------------------------------------------------------------
  -  end function
  +  end function VTK_GEO_XML_UNST_R4_READ
  +
  +
  +  function VTK_GEO_XML_UNST_R8_READ(NN,NC,X,Y,Z,npiece,cf) result(E_IO)
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  !< Function for reading mesh with \b UnstructuredGrid topology (R8P).
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  integer(I4P),           intent(OUT) :: NN       !< number of nodes
  +  integer(I4P),           intent(OUT) :: NC       !< number of cells
  +  real(R8P), allocatable, intent(OUT) :: X(:)     !< x coordinates
  +  real(R8P), allocatable, intent(OUT) :: Y(:)     !< y coordinates
  +  real(R8P), allocatable, intent(OUT) :: Z(:)     !< z coordinates
  +  integer(I4P), optional, intent(IN)  :: npiece   !< Number of the piece to read (by default: 1)
  +  integer(I4P), optional, intent(IN)  :: cf       !< Current file index (for concurrent files IO).
  +  integer(I4P)                        :: E_IO     !< Input/Output inquiring flag: $0$ if IO is done, $> 0$ if IO is not done
  +  integer(I4P)                        :: rf       !< Real file index.
  +  character(len=:),allocatable        :: s_buffer !< Buffer string.
  +  character(len=:), allocatable       :: fmt
  +  character(len=:), allocatable       :: type
  +!  character(len=:), allocatable       :: bindata
  +  integer(I4P)                        :: np, i, offs, N_Byte
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  rf = f
  +  if (present(cf)) then
  +    rf = cf ; f = cf
  +  endif
  +  np = 1_I4P; if (present(npiece)) np = npiece
  +  E_IO = -1_I4P
  +  select case(vtk(rf)%f)         
  +
  +    case(ascii)
  +! not implemented
  +
  +    case(binary)
  +! not implemented
  +
  +    case(raw)
  +
  +      rewind(unit=vtk(rf)%u, iostat=E_IO)
  +      E_IO = move(inside='UnstructuredGrid', to_find='Piece', repeat=np, cf=rf, buffer=s_buffer) ! find the 'np' piece
  +      call get_int(buffer=s_buffer, attrib='NumberOfPoints', val=NN, E_IO=E_IO)
  +      if(E_IO == 0) then
  +        allocate(X(NN), Y(NN), Z(NN), stat=E_IO)
  +        call get_int(buffer=s_buffer, attrib='NumberOfCells', val=NC, E_IO=E_IO)
  +        E_IO = search(inside='Points', to_find='DataArray', with_attribute='Name', of_value='Points', buffer=s_buffer)
  +        call get_int(buffer=s_buffer,  attrib='offset', val=offs, E_IO=E_IO)
  +        call get_char(buffer=s_buffer, attrib='format', val=fmt,  E_IO=E_IO)
  +        call get_char(buffer=s_buffer, attrib='type', val=type,  E_IO=E_IO)
  +        if (trim(adjustlt(Upper_Case(fmt)))/='APPENDED' .or. &
  +            trim(adjustlt(Upper_Case(type)))/='FLOAT64') then
  +          E_IO = -1_I4P 
  +        else
  +          read(unit=vtk(rf)%u, iostat=E_IO, pos = vtk(rf)%ioffset+offs) N_Byte, (X(i), Y(i), Z(i), i=1,NN) !get appended array
  +        endif
  +      endif
  +    end select
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  end function VTK_GEO_XML_UNST_R8_READ
  +
  +  function VTK_GEO_XML_UNST_PACK_R4_READ(NN,NC,XYZ,npiece,cf) result(E_IO)
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  !< Function for reading mesh with \b UnstructuredGrid topology (R8P).
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  integer(I4P),           intent(OUT) :: NN       !< number of nodes
  +  integer(I4P),           intent(OUT) :: NC       !< number of cells
  +  real(R4P), allocatable, intent(OUT) :: XYZ(:,:)   !< Coordinates
  +  integer(I4P), optional, intent(IN)  :: npiece   !< Number of the piece to read (by default: 1)
  +  integer(I4P), optional, intent(IN)  :: cf       !< Current file index (for concurrent files IO).
  +  integer(I4P)                        :: E_IO     !< Input/Output inquiring flag: $0$ if IO is done, $> 0$ if IO is not done
  +  integer(I4P)                        :: rf       !< Real file index.
  +  character(len=:),allocatable        :: s_buffer !< Buffer string.
  +  character(len=:), allocatable       :: fmt
  +  character(len=:), allocatable       :: type
  +!  character(len=:), allocatable       :: bindata
  +  integer(I4P)                        :: np, i, j, offs, N_Byte
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  rf = f
  +  if (present(cf)) then
  +    rf = cf ; f = cf
  +  endif
  +  np = 1_I4P; if (present(npiece)) np = npiece
  +  E_IO = -1_I4P
  +  select case(vtk(rf)%f)         
  +
  +    case(ascii)
  +! not implemented
  +
  +    case(binary)
  +! not implemented
  +
  +    case(raw)
  +
  +      rewind(unit=vtk(rf)%u, iostat=E_IO)
  +      E_IO = move(inside='UnstructuredGrid', to_find='Piece', repeat=np, cf=rf, buffer=s_buffer) ! find the 'np' piece
  +      call get_int(buffer=s_buffer, attrib='NumberOfPoints', val=NN, E_IO=E_IO)
  +      if(E_IO == 0) then
  +        allocate(XYZ(3,NN), stat=E_IO)
  +        call get_int(buffer=s_buffer, attrib='NumberOfCells', val=NC, E_IO=E_IO)
  +        E_IO = search(inside='Points', to_find='DataArray', with_attribute='Name', of_value='Points', buffer=s_buffer)
  +        call get_int(buffer=s_buffer,  attrib='offset', val=offs, E_IO=E_IO)
  +        call get_char(buffer=s_buffer, attrib='format', val=fmt,  E_IO=E_IO)
  +        call get_char(buffer=s_buffer, attrib='type', val=type,  E_IO=E_IO)
  +        if (trim(adjustlt(Upper_Case(fmt)))/='APPENDED' .or. &
  +            trim(adjustlt(Upper_Case(type)))/='FLOAT32') then
  +          E_IO = -1_I4P 
  +        else
  +          read(unit=vtk(rf)%u, iostat=E_IO, pos = vtk(rf)%ioffset+offs) N_Byte, ((XYZ(i,j),i=1,3),j=1,NN) !get appended array
  +        endif
  +      endif
  +    end select
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  end function VTK_GEO_XML_UNST_PACK_R4_READ
  +
  +
  +  function VTK_GEO_XML_UNST_PACK_R8_READ(NN,NC,XYZ,npiece,cf) result(E_IO)
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  !< Function for reading mesh with \b UnstructuredGrid topology (R8P).
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  integer(I4P),           intent(OUT) :: NN       !< number of nodes
  +  integer(I4P),           intent(OUT) :: NC       !< number of cells
  +  real(R8P), allocatable, intent(OUT) :: XYZ(:,:)   !< Coordinates
  +  integer(I4P), optional, intent(IN)  :: npiece   !< Number of the piece to read (by default: 1)
  +  integer(I4P), optional, intent(IN)  :: cf       !< Current file index (for concurrent files IO).
  +  integer(I4P)                        :: E_IO     !< Input/Output inquiring flag: $0$ if IO is done, $> 0$ if IO is not done
  +  integer(I4P)                        :: rf       !< Real file index.
  +  character(len=:),allocatable        :: s_buffer !< Buffer string.
  +  character(len=:), allocatable       :: fmt
  +  character(len=:), allocatable       :: type
  +!  character(len=:), allocatable       :: bindata
  +  integer(I4P)                        :: np, i, j, offs, N_Byte
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  rf = f
  +  if (present(cf)) then
  +    rf = cf ; f = cf
  +  endif
  +  np = 1_I4P; if (present(npiece)) np = npiece
  +  E_IO = -1_I4P
  +  select case(vtk(rf)%f)         
  +
  +    case(ascii)
  +! not implemented
  +
  +    case(binary)
  +! not implemented
  +
  +    case(raw)
  +
  +      rewind(unit=vtk(rf)%u, iostat=E_IO)
  +      E_IO = move(inside='UnstructuredGrid', to_find='Piece', repeat=np, cf=rf, buffer=s_buffer) ! find the 'np' piece
  +      call get_int(buffer=s_buffer, attrib='NumberOfPoints', val=NN, E_IO=E_IO)
  +      if(E_IO == 0) then
  +        allocate(XYZ(3,NN), stat=E_IO)
  +        call get_int(buffer=s_buffer, attrib='NumberOfCells', val=NC, E_IO=E_IO)
  +        E_IO = search(inside='Points', to_find='DataArray', with_attribute='Name', of_value='Points', buffer=s_buffer)
  +        call get_int(buffer=s_buffer,  attrib='offset', val=offs, E_IO=E_IO)
  +        call get_char(buffer=s_buffer, attrib='format', val=fmt,  E_IO=E_IO)
  +        call get_char(buffer=s_buffer, attrib='type', val=type,  E_IO=E_IO)
  +        if (trim(adjustlt(Upper_Case(fmt)))/='APPENDED' .or. &
  +            trim(adjustlt(Upper_Case(type)))/='FLOAT64') then
  +          E_IO = -1_I4P 
  +        else
  +          read(unit=vtk(rf)%u, iostat=E_IO, pos = vtk(rf)%ioffset+offs) N_Byte, ((XYZ(i,j),i=1,3),j=1,NN) !get appended array
  +        endif
  +      endif
  +    end select
  +  !---------------------------------------------------------------------------------------------------------------------------------
  +  end function VTK_GEO_XML_UNST_PACK_R8_READ
   
   
     function VTK_GEO_XML_STRG_1DA_R4_READ(nx1,nx2,ny1,ny2,nz1,nz2,NN,X,Y,Z,npiece,cf) result(E_IO)
  
positive.code:       endif
      end select
    !---------------------------------------------------------------------------------------------------------------------------------
    end function
  
  
    function VTK_GEO_XML_RECT_R4_READ(nx1,nx2,ny1,ny2,nz1,nz2,X,Y,Z,npiece,cf) result(E_IO)
  
negative.code:       E_IO = move(inside='RectilinearGrid', to_find='Piece', repeat=np, cf=rf, buffer=s_buffer) ! find the 'np' piece
        call get_char(buffer=s_buffer, attrib='Extent', val=aux, E_IO=E_IO)
        if(E_IO == 0) then      
          read(aux,*, iostat=E_IO) nx1,nx2,ny1,ny2,nz1,nz2
          if(E_IO == 0) then      
            NNp = (nx2-nx1+1)*(ny2-ny1+1)*(nz2-nz1+1)
            allocate(X(NNp), Y(NNp), Z(NNp), stat=E_IO)        
            E_IO = search(inside='Coordinates', to_find='DataArray', with_attribute='Name', of_value='X', buffer=s_buffer)
  