C 
      SUBROUTINE SOFT
     *           (NCLASS,NA1,NA2,NA3,NB1,NB2,NB3,WEIGHT,NM,INDM,CS,MB,B)
C
      INTEGER NCLASS,NA1,NA2,NA3,NB1,NB2,NB3,NM,INDM(*),MB
      REAL    WEIGHT,CS(*),B(MB)
C
C This subroutine accumulates the prior subjective information
C covariance matrix describing the smoothness of the functions
C interpolated by means of subroutines of the file 'val.for'.
C
C Input:
C     NCLASS..Number of classes required.
C     NA1,NA2,NA3... Orders of X1,X2,X3-partial derivatives of the first
C             basis function in the product being integrated.
C     NB1,NB2,NB3... Orders of X1,X2,X3-partial derivatives of the
C             second basis function in the product being integrated.
C     WEIGHT..Weighting factor corresponding to the derivatives.
C     CS...   Symmetric matrix to be regularized:
C             Contribution corresponding to the derivatives will be
C             added to the matrix CS if WEIGHT.NE.0.
C             Not used if WEIGHT=0.
C     MB...   Dimension of working array B.
C     B...    Working array.  Auxiliary storage locations for partial
C             covariance matrices.  For the comments on its minimum
C             dimension refer to the destription of error
C             372.  Not used if WEIGHT=0.
C
C Output:
C     NM...   Number of the unknown model parameters.
C     INDM... Indices of the unknown model parameters.  Integer array of
C             dimension NM.
C     CS...   Resulting symmetric subjective prior information
C             covariance matrix (input matrix increased by WEIGHT times
C             the product of the derivatives averaged over the B-spline
C             grid).  Real array of dimension NM*(NM+1)/2.
C             Unchanged input if WEIGHT=0.
C
C Common block:
      INCLUDE 'val.inc'
C     val.inc
C None of the storage locations of the common block are altered.
C
C Subroutines and external functions required:
      EXTERNAL SPSP
C
C Date: 1997, September 30
C Coded by Ludek Klimes
C
C-----------------------------------------------------------------------
C
      INTEGER ICLASS,JGROUP,LFUNCT,MFUNCT,JFUNCT,LADR,MADR,IADR,NVAR
      INTEGER NX(3),NX1,NX2,NX3
      EQUIVALENCE (NX(1),NX1),(NX(2),NX2),(NX(3),NX3)
      INTEGER JADR(7),JADR1,JADR2,JADR3,JADR4,JADR5,JADR6,JADR7
      EQUIVALENCE (JADR(1),JADR1),(JADR(2),JADR2),(JADR(3),JADR3)
      EQUIVALENCE (JADR(4),JADR4),(JADR(5),JADR5),(JADR(6),JADR6)
      EQUIVALENCE (JADR(7),JADR7)
      INTEGER I1,I2,I3,II,I,J,N
      REAL WDENS,SIGMA
C
C     B...    Auxiliary storage locations for partial covariance
C             matrices.
C     ICLASS..Index of the class.
C     JGROUP..Address of a current group.
C     LFUNCT,MFUNCT,JFUNCT... Addresses of the first, last and arbitrary
C             functions of the group.
C     LADR,MADR,IADR... Addresses of the first, last and arbitrary
C             parameters of the current function.
C     NVAR... Number of the independent variables A1, A2, A3 of the
C             interpolated function W.
C     NX=(NX1,NX2,NX3)... Numbers of grid lines.
C     JADR=(JADR1,JADR2,JADR3,JADR4,JADR5,JADR6,JADR7)... Addresses of
C             parameters describing the interpolated function (grid
C             coordinates, B-spline coefficients, B-spline basis
C             functions).
C     I1,I2,I3,I,J,N... Local auxiliary variables.
C     WDENS...Weighting factor divided by volume (area,length).
C     SIGMA...Tension factor.
C
C.......................................................................
C
      NM=0
      IF(NCLASS.LT.1.OR.IPAR(0).LT.NCLASS) THEN
C       371
        CALL ERROR('371 in SOFT: Incorrect number of classes')
C       The number of classes required is zero, negative or greater than
C       the number of classes defined.
      END IF
C     Loop for classes:
      DO 92 ICLASS=1,NCLASS
C       Loop for groups of the current class:
        DO 91 JGROUP=IPAR(ICLASS-1)+1,IPAR(ICLASS)
          LFUNCT=IPAR(JGROUP-1)+1
          MFUNCT=IPAR(JGROUP)
          MADR  =IPAR(LFUNCT-1)
C
C         Loop for functions of the current group:
          DO 90 JFUNCT=LFUNCT,MFUNCT
C           Starting and end addresses of the parameters describing the
C           function
            LADR=MADR+1
            MADR=IPAR(JFUNCT)
            IF(LADR.LE.MADR) THEN
              WDENS=WEIGHT
C             Tension factor
              SIGMA=RPAR(LADR+5)
C
C             The number, types and values of the independent variables
C             of the interpolated function:
C             Initial address
              IADR=LADR+6
C             Initial number of the independent variables
              NVAR=0
              NX1=1
              NX2=1
              NX3=1
              JADR1=0
              JADR2=0
              JADR3=0
              JADR4=0
C             Loop for the possible independent variables:
              DO 20 I=LADR+2,LADR+4
C               Type of the possible independent variable:
                J=IPAR(I)
                IF(J.GT.0) THEN
                  N=IABS(IPAR(IADR))
                  IF(N.GE.2) THEN
                    NVAR=NVAR+1
                    NX(NVAR)=N
                  ELSE IF(N.EQ.1) THEN
                    JADR(NVAR+1)=JADR(NVAR+1)+1
                  END IF
                  IADR=IADR+1
                END IF
   20         CONTINUE
              JADR5=0
              JADR6=0
              JADR7=0
C
C             Interpolated function W:
              JADR1=IADR+JADR1
              IF(NVAR.LE.0) THEN
C               No independent variable.
                JADR4=JADR1
              ELSE
                JADR2=JADR1+NX1+JADR2
                WDENS=WDENS/(RPAR(JADR1+NX1-1)-RPAR(JADR1))
                IF(NVAR.EQ.1) THEN
C                 One independent variable:
                  JADR4=JADR2
                  JADR5=JADR2+NX1
                ELSE
                  JADR3=JADR2+NX2+JADR3
                  WDENS=WDENS/(RPAR(JADR2+NX2-1)-RPAR(JADR2))
                  IF(NVAR.EQ.2) THEN
C                   Two independent variables:
                    JADR4=JADR3
                    JADR5=JADR3+NX1*NX2
                    JADR6=JADR4+5*NX1
                  ELSE
C                   Three independent variables:
                    JADR4=JADR3+NX3+JADR4
                    JADR5=JADR4+NX1*NX2*NX3
                    JADR6=JADR5+5*NX1
                    JADR7=JADR6+5*NX2
                    WDENS=WDENS/(RPAR(JADR3+NX3-1)-RPAR(JADR3))
                  END IF
                END IF
              END IF
C             WDENS is weighting factor divided by volume (area,length).
C
              JADR4=JADR4-1
              N=NX1*NX2*NX3
              IF(WEIGHT.EQ.0.) THEN
                DO 72 I2=1,N
                  INDM(NM+I2)=JADR4+I2
   72           CONTINUE
                NM=NM+N
              ELSE
                IF(NA1.EQ.NB1.AND.NA2.EQ.NB2.AND.NA3.EQ.NB3) THEN
                  I=N*(N+1)/2
                  J=0
                ELSE
                  I=N*N
                  J=1
                END IF
                IF(NA1.EQ.NB1) THEN
                  I1=NX1*(NX1+1)/2
                ELSE
                  I1=NX1*NX1
                END IF
                IF(NA2.EQ.NB2) THEN
                  I2=NX2*(NX2+1)/2
                ELSE
                  I2=NX2*NX2
                END IF
                IF(NA3.EQ.NB3) THEN
                  I3=NX3*(NX3+1)/2
                ELSE
                  I3=NX3*NX3
                END IF
                I1=I1+I
                I2=I2+I1
                I3=I3+I2
                IF(I3.GT.MB) THEN
C                 372
                  CALL ERROR('372 in SOFT: Insufficient working memory')
C                 The dimension MB of the buffer B is not sufficient,
C                 see the above command lines, where NX1,NX2,NX3 are
C                 the dimensions of the grid for spline interpolation.
C                 MB and B(MB) are the dummy arguments of this
C                 subroutine.  If the actual argument is allocated in
C                 array RAM(MRAM) declared in include file
C                 ram.inc,
C                 you may wish to increase MRAM.
                END IF
                DO 80 II=1,I
                  B(II)=0.
   80           CONTINUE
                CALL SPSP(NA1,NA2,NA3,NB1,NB2,NB3,NX1,NX2,NX3,
     *                    RPAR(JADR1),RPAR(JADR2),RPAR(JADR3),
     *                    RPAR(JADR5),RPAR(JADR6),RPAR(JADR7),
     *                    SIGMA,WDENS,B(1),B(I+1),B(I1+1),B(I2+1))
                I=NM*(NM+1)/2
                DO 82 I2=1,N
                  INDM(NM+I2)=JADR4+I2
                  I=I+NM
                  DO 81 I1=1,I2
                    I=I+1
                    IF(J.EQ.0) THEN
                      CS(I)=CS(I)+B(I2*(I2-1)/2+I1)
                    ELSE
                      CS(I)=CS(I)+B(N*(I2-1)+I1)+B(N*(I1-1)+I2)
                    END IF
   81             CONTINUE
   82           CONTINUE
                NM=NM+N
C               Contribution corresponding to the derivatives
C               is added to the matrix CS.
              END IF
C
            END IF
   90     CONTINUE
   91   CONTINUE
   92 CONTINUE
C     End of loops for functions.
C
      RETURN
      END
C
C=======================================================================
C