#! /bin/sh
#
#    Last change: Mihai Dima 2001
#
if [ $# -lt 3 ]
then
  echo " "
  echo "            *****     NRMSTDV     *****      "
  echo " "
  echo "The program calculates the standard deviation field (file3)"
  echo "and normalizes the input data (file1) in output file (file2)"
  echo " "
  echo "                     Syntax: nrmstdv <file1> <file2> <file3>"
  echo " "
  exit 0
else
 echo " "
 echo "            *****     NRMSTDV     *****      "
 echo " "
 echo "file1="$1
 echo "file2="$2
 echo "file3="$3
 cp $1 fort.1
fi
cat > nrmstdv.f << EOF
      PROGRAM NRMSTDV
      PARAMETER(UNDEF=0.9E+10)
      INTEGER H(4),I,J

      REAL F(70000),SX(70000),SMED(70000)
      REAL STDV(70000),SN(70000)


      NCHAN=1
      CALL PARLES (NCHAN,IDIM,NSET)
     
      READ (1) H
      IDIM=H(4)      
      
      DO I=1,IDIM
        SMED(I)=0.0
        SX(I)=0.0
        SN(I)=0.0 
      ENDDO

      REWIND 1
    

60    READ (1,END=61) H
      READ (1) (F(J),J=1,IDIM)
      DO I=1,H(4)
        IF (F(I).NE.UNDEF) THEN
          SMED(I)=SMED(I)+F(I)
          SN(I)=SN(I)+1.0
        ENDIF
      ENDDO
      GOTO 60

61    CONTINUE

      DO I=1,H(4)
         IF (SN(I).NE.0.0) THEN
            SMED(I)=SMED(I)/SN(I)   
         ENDIF
      ENDDO

      REWIND 1
      

50    READ (1,END=51) H
      READ (1) (F(J),J=1,IDIM)
      DO 43 I=1,H(4)
        IF (SMED(I).NE.0.0) THEN
          SX(I)=SX(I)+(F(I)-SMED(I))*(F(I)-SMED(I))
        ELSE
          SX(I)=UNDEF  
        ENDIF
43    CONTINUE
      GOTO 50

51    CONTINUE

      DO I=1,H(4)
         IF (SX(I).NE.UNDEF) THEN
            STDV(I)=SQRT(SX(I)/(SN(I)-1.0))
         ELSE
            STDV(I)=UNDEF      
         ENDIF
      ENDDO


      WRITE(3) H(1),H(2),H(3),IDIM
      WRITE(3) (STDV(I),I=1,IDIM)


      REWIND 1

70    READ (1,END=71) H
      IDIM=H(4)
      READ (1) (F(J),J=1,IDIM)
      DO I=1,IDIM
        IF(F(I).NE.UNDEF) THEN
          F(I)=F(I)/STDV(I)
        ENDIF
      ENDDO

      WRITE(2) H
      WRITE(2) (F(J),J=1,IDIM)
      GOTO 70

71    CONTINUE
        

      CLOSE(3)
      CLOSE(2)
      CLOSE(1)

      PRINT *,'REGULAR END!'

      END
C******************************************************************
      SUBROUTINE PARLES (NCHAN,NDIM,NSET)
      IMPLICIT REAL (A-H,O-Z)

      INTEGER I1,I2,I3

      NSET=0

      REWIND NCHAN
      READ(NCHAN,END=10) ITIME,INAME,ILEV,ISIZE
      GOTO 20
   10 PRINT *,'INPUT FILE ON TAPE',NCHAN,' IS EMPTY'
      STOP 'MISTAKE'
   20 REWIND NCHAN

      NDIM=0  
    1 READ (NCHAN,END=100) I1,I2,I3,I4
      NDIM=MAX0(NDIM,I4)
      READ (NCHAN,END=101)
      NSET=NSET+1
      GOTO1
      
  101 PRINT *,' MORE HEADER THEN DATA '     
  100 PRINT 1000,NCHAN,NDIM,NSET
      REWIND NCHAN
      RETURN    
 1000 FORMAT (1X,' TAPE=',I2,',  NDIM=',I5,',  NSET=',I5)
      END
EOF
f77 nrmstdv.f -o nrmstdv.x
nrmstdv.x $1 $2 $3
cp fort.2 $2
cp fort.3 $3
rm fort.1 fort.2 fort.3
rm nrmstdv.x nrmstdv.f
exit