next up previous contents index
Next: Initialization file Up: A Template Task Previous: A Template Task

Source code

        PROGRAM COMBINE
C----------------------------------------------------------------------
C GILDAS Combine in different ways two input images
C        (or data cubes)...
C----------------------------------------------------------------------
        USE IMAGE_DEF
        INCLUDE 'inc:format.inc'
        INCLUDE 'inc:errcod.inc'
*
        CHARACTER*80 NAMEX,NAMEY,NAMEZ, CODE*20
        LOGICAL ERROR
        REAL AY,AZ,TY,TZ,B,C
        REAL, ALLOCATABLE :: DX(:,:), DY(:,:), DZ(:)
        INTEGER I, J, N, M, IER, LENC
*
        CALL GILDAS_OPEN
        CALL GILDAS_CHAR('Z_NAME$',NAMEZ)
        CALL GILDAS_REAL('Z_FACTOR$',AZ,1) 
        CALL GILDAS_REAL('Z_MIN$',TZ,1)
        CALL GILDAS_CHAR('Y_NAME$',NAMEY)
        CALL GILDAS_REAL('Y_FACTOR$',AY,1) 
        CALL GILDAS_REAL('Y_MIN$',TY,1)
        CALL GILDAS_CHAR('X_NAME$',NAMEX)
        CALL GILDAS_REAL('BLANKING$',B,1)
        CALL GILDAS_REAL('OFFSET$',C,1)
        CALL GILDAS_CHAR('FUNCTION$',CODE) 
        CALL GILDAS_CLOSE
*
        N = LENC(NAMEZ)
        IF (N.EQ.0) GOTO 100
        CALL GILDAS_NULL(HZ)
        CALL SIC_PARSEF(NAMEZ(1:N),HZ%FILE,' ','.gdf')
        CALL GDF_READ_HEADER(HZ,ERROR)
        IF (ERROR) THEN
          WRITE(6,*) 'F-COMBINE,  Cannot read input file'
          GOTO 100
        ENDIF
*
        N = LENC(NAMEY)
        IF (N.EQ.0) GOTO 100
        CALL GILDAS_NULL(HY)
        CALL SIC_PARSEF(NAMEY(1:N),HY%FILE,' ','.gdf')
        CALL GDF_READ_HEADER(HY,ERROR)
        IF (ERROR) THEN
          WRITE(6,*) 'F-COMBINE,  Cannot read input file'
          GOTO 100
        ENDIF
*
        IF (HZ%GIL%EVAL.GE.0.0) HZ%GIL%EVAL = 
      & 	MAX(HZ%GIL%EVAL,ABS(HZ%GIL%BVAL*1E-7))
        IF (HY%GIL%EVAL.GE.0.0) HY%GIL%EVAL = 
      & 	MAX(HY%GIL%EVAL,ABS(HY%GIL%BVAL*1E-7))
*
* Check input dimensions
        DO I=1,4
          IF (HY%GIL%DIM(I).NE.HZ%GIL%DIM(I)) THEN
            N = 1
            DO J=I,4
              N = N*HZ%GIL%DIM(J)
            ENDDO
            IF (N.NE.1) THEN
              WRITE(6,*) 'F-COMBINE,  Input images are non coincident'
              GOTO 100
            ELSE
              WRITE(6,*) 'W-COMBINE,  Combining a cube with a plane'
            ENDIF
          ENDIF
        ENDDO
*
        CALL GDF_COPY_HEADER(HY,HX)
        N = LENC(NAMEX)
        IF (N.EQ.0) GOTO 100
        CALL SIC_PARSEF(NAMEX(1:N),HX%FILE,' ','.gdf')
        HX%GIL%BLAN = 2
        HX%GIL%BVAL = B
        HX%GIL%EVAL = 0.0
        HX%GIL%EXTR = 0
*
* Allocate the arrays. Note that the allocated arrays do not conform
* to the shape of the images: DZ is allocated as a 1-D array, DX,DY
* as 2-D arrays, possibly of second dimension 1.
*
        N = HZ%LOCA%SIZE
        M = HX%LOCA%SIZE/HZ%LOCA%SIZE
        ALLOCATE(DX(N,M),DY(N,M),DX(N),STAT=IER)
        IF (IER.NE.0) THEN
          WRITE(6,*) 'F-COMBINE,  Input images are non coincident'
          GOTO 100
        ENDIF
*
* Read the input data
        CALL GDF_READ_DATA(HZ,DZ,ERROR)
        CALL GDF_READ_DATA(HY,DY,ERROR)
*
        IF (CODE.EQ.'ADD') THEN
          CALL ADD002(DZ,DY,DX,
     +      N,M,
     +      HZ%GIL%BVAL,HZ%GIL%EVAL,AZ,TZ,
     +      HY%GIL%BVAL,HY%GIL%EVAL,AY,TY,
     +      HX%GIL%BVAL,C)
        ELSEIF (CODE.EQ.'DIVIDE') THEN
          CALL DIV002(DZ,DY,DX,
     +      N,M,
     +      HZ%GIL%BVAL,HZ%GIL%EVAL,AZ,TZ,
     +      HY%GIL%BVAL,HY%GIL%EVAL,AY,TY,
     +      HX%GIL%BVAL,C)
        ELSEIF (CODE.EQ.'MULTIPLY') THEN
          CALL MUL002(DZ,DY,DX,
     +      N,M,
     +      HZ%GIL%BVAL,HZ%GIL%EVAL,AZ,TZ,
     +      HY%GIL%BVAL,HY%GIL%EVAL,AY,TY,
     +      HX%GIL%BVAL,C)
        ELSEIF (CODE.EQ.'OPTICAL_DEPTH') THEN
          CALL OPT002(DZ,DY,DX,
     +      N,M,
     +      HZ%GIL%BVAL,HZ%GIL%EVAL,AZ,TZ,
     +      HY%GIL%BVAL,HY%GIL%EVAL,AY,TY,
     +      HX%GIL%BVAL,C)
        ELSE
            WRITE(6,*) 'Invalid operation code '//CODE
            GOTO 100
        ENDIF
*
* Write ouput file 
        CALL GDF_WRITE_IMAGE(HX,DX,ERROR)
*
        STOP 'S-COMBINE,  Successful completion'
*
98      WRITE(6,*) 'F-COMBINE,  Missing parameter'
100     CALL SYSEXI (FATALE)
        END
*
**
        SUBROUTINE ADD002(Z,Y,X,N,M,BZ,EZ,AZ,TZ,BY,EY,AY,TY,BX,C)
C----------------------------------------------------------------------
C TASK  Internal routine
C       Linear combination of input arrays
C       X = Ay*Y + Az*Z + C
C Arguments
C       Z               R*4(*)  Input array (N)
C       Y               R*4(*)  Input array (N,M)
C       X               R*4(*)  Output array (N,M)
C       N,M             I       Dimensions of arrays
C       BX,BY,BZ        R*4     Blanking values
C       EY,EZ           R*4     Tolerance on blanking
C       AZ,AY           R*4     Multiplicative factor of array Z, Y
C       TZ,TY           R*4     Threshold on Z,Y
C       C               R*4     Additive constant 
C----------------------------------------------------------------------
        INTEGER N,M
        REAL Z(N),Y(N,M),X(N,M),BX,BY,BZ,EY,EZ,AZ,AY,TY,TZ,C
        INTEGER I,K
*
        DO K=1,M
          DO I=1,N
            IF (ABS(Z(I)-BZ).GT.EZ .AND. ABS(Y(I,K)-BY).GT.EY 
     +      .AND. Z(I).GT.TZ .AND. Y(I,K).GT.TY) THEN
               X(I,K) = AY*Y(I,K) + AZ*Z(I)  +  C 
            ELSE
               X(I,K) = BX
            ENDIF
          ENDDO
        ENDDO
        END



Gildas manager
2002-06-07