decompress.f90

Go to the documentation of this file.

!-----------------------------------------------------------------
!--------------- special set of characters for RCS information
!-----------------------------------------------------------------
! $Source: /home/cvsroot/mesonh/libtools/lib/COMPRESS/src/decompress.f90,v $ $Revision: 1.1 $ $Date: 2005/04/12 15:06:20 $
!-----------------------------------------------------------------
SUBROUTINE get_compheader(KTAB,SIZEKTAB,KNBELT,KTYPECOD)

INTEGER, INTENT(IN) :: SIZEKTAB
INTEGER(KIND=8), DIMENSION(SIZEKTAB), INTENT(IN) :: KTAB
INTEGER, INTENT(OUT) :: KNBELT    ! size of decompressed array
INTEGER, INTENT(OUT) :: KTYPECOD  ! code for compression type

 CHARACTER(LEN=8) :: STRKEY

INTEGER :: INTCHAR
INTEGER :: JI

 CALL set_extractidx(0,0)
! extract string header 
DO ji=1,8
  CALL extract_bbuff(ktab,8,intchar)
  strkey(ji:ji) = char(intchar)
END DO

! Treat array if it is compressed
IF (strkey == 'COMPRESS') THEN
  CALL extract_bbuff(ktab,32,ktypecod)
  CALL extract_bbuff(ktab,32,knbelt)
ELSE
  knbelt    =-1
  ktypecod = 0
END IF

END SUBROUTINE get_compheader

SUBROUTINE decompress_field(XTAB,NBELT,COMPTAB,NBCOMPELT,CODINGTYPE)
USE modd_comppar
USE mode_searchgrp

IMPLICIT NONE 
INTEGER,                                INTENT(IN)  :: NBELT 
INTEGER,                                INTENT(IN)  :: NBCOMPELT 
REAL   (KIND=8),DIMENSION(NBELT),TARGET,INTENT(OUT) :: XTAB
INTEGER(KIND=8),DIMENSION(NBCOMPELT),   INTENT(IN)  :: COMPTAB
INTEGER,                                INTENT(IN)  :: CODINGTYPE

INTEGER,DIMENSION(:), ALLOCATABLE  :: ITAB
LOGICAL,DIMENSION(:), ALLOCATABLE  :: GMASK

REAL :: XREF, XCOEFF
INTEGER :: INBLEV
INTEGER :: ILEVNBELT
INTEGER :: JI
INTEGER :: IND1, IND2
INTEGER :: IDIMX,IDIMY
INTEGER :: IEXTCOD
REAL(KIND=8),DIMENSION(:),POINTER  :: XPTRTAB
REAL :: XMIN,XMAX

SELECT CASE (codingtype)
CASE (jpcstencod)
  CALL extract_bbuff(comptab,32,xref)
  xtab(:) = xref

CASE (jpsopencod)
  CALL extract_bbuff(comptab,32,idimx)
  CALL extract_bbuff(comptab,32,idimy)
  ilevnbelt = idimx * idimy
  inblev = nbelt/(ilevnbelt)
  ALLOCATE(itab(ilevnbelt))
  DO ji=1,inblev
    ind1=(ji-1)*ilevnbelt+1
    ind2=ji*ilevnbelt
    xptrtab=>xtab(ind1:ind2)
    IF (lpdebug) print *,'######   Decompress(SOPENCOD) LEVEL ',ji,'######'
    CALL extract_bbuff(comptab,32,xref)
    CALL extract_bbuff(comptab,32,xcoeff)
    CALL extractintarray(itab)
    CALL decomp_fop(xptrtab,itab,xref,xcoeff)
  END DO
  
CASE (jpextencod)
  CALL extract_bbuff(comptab,32,idimx)
  CALL extract_bbuff(comptab,32,idimy)
  ilevnbelt = idimx * idimy
  inblev = nbelt/(ilevnbelt)
  ALLOCATE(itab(ilevnbelt))
  ALLOCATE(gmask(ilevnbelt))
  DO ji=1,inblev

    IF (lpdebug) print *,'###### Decompress(EXTENCOD) LEVEL ',ji,'######'
    ind1=(ji-1)*ilevnbelt+1
    ind2=ji*ilevnbelt
    xptrtab=>xtab(ind1:ind2)
    !
    CALL extract_bbuff(comptab,3,iextcod)
    IF (iextcod == jpother) THEN
      CALL extract_bbuff(comptab,3,iextcod)
      iextcod = iextcod + 8
    END IF
    IF (lpdebug) print *, "IEXTCOD = ",iextcod
    SELECT CASE(iextcod)
    CASE(jplog)
      ! Conversion to log values of original data 0<=x<1
      CALL extract_bbuff(comptab,32,xref)
      CALL extract_bbuff(comptab,32,xcoeff)
      CALL extractintarray(itab)
      gmask(:) = .true.
      WHERE (itab == 0)
        gmask = .false.
        xptrtab = 0.0
      END WHERE
      CALL decomp_fop(xptrtab,itab,xref,xcoeff,gmask,1)
      WHERE(gmask)
        xptrtab = exp(xptrtab)
      END WHERE
      
    CASE(jpconst)
      ! constant value array
      CALL extract_bbuff(comptab,32,xref)
      xptrtab(:) = xref
      IF (lpdebug) print *,"  CONST value=",xref

    CASE(jp2val)
      ! 2 different values in array
      CALL extract_bbuff(comptab,32,xmin)
      CALL extract_bbuff(comptab,32,xmax)
      CALL extractintarray(itab)
      WHERE (itab == 0)
        xptrtab = xmin
      ELSEWHERE
        xptrtab = xmax
      END WHERE
      IF (lpdebug) print *,"  2 values:",xmin,xmax
      
    CASE(jp3val)
      ! 3 different values in array
      CALL extract_bbuff(comptab,32,xmin)
      CALL extract_bbuff(comptab,32,xref)
      CALL extract_bbuff(comptab,32,xmax)
      CALL extractintarray(itab)
      WHERE (itab == 0)
        xptrtab = xmin
      ELSEWHERE
        xptrtab = xref
      END WHERE
      WHERE (itab == 2) xptrtab = xmax
      IF (lpdebug) print *,"  3 values:",xmin,xref,xmax

    CASE(jpnorm)
      ! same as JPSOPENCOD
      CALL extract_bbuff(comptab,32,xref)
      CALL extract_bbuff(comptab,32,xcoeff)
      CALL extractintarray(itab)
      CALL decomp_fop(xptrtab,itab,xref,xcoeff)
      IF (lpdebug) print *,"  normal, XREF/XCOEFF = ",xref,xcoeff 

    CASE(jpminexcl)
      ! Min value is isolated
      CALL extract_bbuff(comptab,32,xmin)
      CALL extract_bbuff(comptab,32,xref)
      CALL extract_bbuff(comptab,32,xcoeff)
      CALL extractintarray(itab)
      gmask(:) = .true.
      WHERE (itab == 0)
        gmask = .false.
        xptrtab = xmin
      END WHERE
      CALL decomp_fop(xptrtab,itab,xref,xcoeff,gmask,1)
      IF (lpdebug) print *,"  Min exclus, MIN/XREF/XCOEFF = ",xmin,xref,xcoeff

    CASE(jpmaxexcl)
      ! Max value is isolated
      CALL extract_bbuff(comptab,32,xmax)
      CALL extract_bbuff(comptab,32,xref)
      CALL extract_bbuff(comptab,32,xcoeff)
      CALL extractintarray(itab)
      gmask(:) = .true.
      WHERE (itab == 65535)
        gmask = .false.
        xptrtab = xmax
      END WHERE
      CALL decomp_fop(xptrtab,itab,xref,xcoeff,gmask,0)   
      IF (lpdebug) print *,"  Max exclus, MAX/XREF/XCOEFF = ",xmax,xref,xcoeff

    CASE(jpminmaxexcl)
      ! Min&Max value are isolated
      CALL extract_bbuff(comptab,32,xmin)        
      CALL extract_bbuff(comptab,32,xmax)
      CALL extract_bbuff(comptab,32,xref)
      CALL extract_bbuff(comptab,32,xcoeff)
      CALL extractintarray(itab)
      gmask(:) = .true.
      WHERE (itab == 0)
        gmask = .false.
        xptrtab = xmin
      END WHERE
      WHERE (itab == 65535)
        gmask = .false.
        xptrtab = xmax
      END WHERE
      CALL decomp_fop(xptrtab,itab,xref,xcoeff,gmask,1)
      IF (lpdebug) print *,"  Min et Max exclus, MIN/MAX/XREF/XCOEFF = ",&
           &xmin,xmax,xref,xcoeff
    END SELECT
  END DO
  
CASE DEFAULT
  print *,'Error in CODINGTYPE : program aborted'
  stop
END SELECT

CONTAINS 

SUBROUTINE decomp_fop(PTAB,KTAB,PREF,PCOEFF,OMASK,KINDCOR)
REAL(KIND=8), DIMENSION(:), INTENT(INOUT) :: PTAB 
! Attention: avec le compilateur PGF, utiliser INTENT(OUT) provoque une recopie
! complete du tableau dans PTAB (avec ecrasement possible des valeurs 
! presentes a l'appel de la procedure). Le phenomene est genant lorsque
! DECOMP_FOP ne calcule que sur une portion de PTAB (valeurs min et/ou max 
! sont presentes). En declarant PTAB en INOUT, les valeurs en entree de la routine
! sont conservees si elles n'ont pas ete modifiees.

INTEGER,      DIMENSION(:), INTENT(IN) :: KTAB 
REAL, INTENT(IN) :: PREF
REAL, INTENT(IN) :: PCOEFF
LOGICAL, DIMENSION(:),INTENT(IN),OPTIONAL :: OMASK
INTEGER,INTENT(IN),OPTIONAL  :: KINDCOR ! 1 if Min value is isolated, 0 otherwise

INTEGER :: INDCOR

IF (.NOT. PRESENT(kindcor)) THEN
  indcor = 0
ELSE
  indcor = kindcor
END IF
  
IF (PRESENT(omask)) THEN
  WHERE (omask)
    ptab(:) = pcoeff*(ktab(:)-indcor)+pref
  END WHERE
ELSE
  IF (pcoeff == 0.0) THEN
    ptab(:) = pref
  ELSE
    ptab(:) = pcoeff*ktab(:)+pref
  END IF
END IF

END SUBROUTINE decomp_fop

SUBROUTINE extractintarray(KTAB)
INTEGER,DIMENSION(:),INTENT(OUT) :: KTAB
!
! COMPTAB, IDIMX and IDIMY  are defined in the calling routine
!
INTEGER :: NBGRP
INTEGER :: IBE
INTEGER :: CPT
INTEGER :: JJ
INTEGER :: ALONE
INTEGER :: NBITCOD,IMIN
INTEGER :: GELT
INTEGER :: JELT
INTEGER :: IEPS

 CALL extract_bbuff(comptab,32,nbgrp)
!      PRINT *,'Nbre de groupes =',NBGRP
 CALL extract_bbuff(comptab,5,ibe)
!      PRINT *,'Nbre de bits pour coder le nombre d''elements:',IBE
cpt = 1
DO jj=1,nbgrp
  !      PRINT *,'Groupe ',JJ,' : '
  CALL extract_bbuff(comptab,1,alone)
  CALL extract_bbuff(comptab,16,imin)
  !      PRINT *,'IREF=',IMIN
  
  IF (alone == 1) THEN
    ! 1 seul elt dans le groupe
    !        PRINT *,'--> un seul element dans le groupe'
    ktab(cpt)=imin
    cpt=cpt+1
  ELSE
    CALL extract_bbuff(comptab,4,nbitcod)
    CALL extract_bbuff(comptab,ibe,gelt)
    !        PRINT *,'--> ',GELT,' elts, codage ecart sur ',nbitcod,'bits'
    IF (nbitcod > 0) THEN
      DO jelt=1,gelt
        CALL extract_bbuff(comptab,nbitcod,ieps)
        ktab(cpt) = imin+ieps
        cpt=cpt+1
      END DO
    ELSE
      ktab(cpt:cpt+gelt-1) = imin
      cpt = cpt+gelt
    END IF
  END IF
END DO
 CALL invertcol(ktab,idimx,idimy)        
END SUBROUTINE extractintarray

END SUBROUTINE decompress_field