home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Frostbyte's 1980s DOS Shareware Collection
/
floppyshareware.zip
/
floppyshareware
/
DOOG
/
PCSSP1.ZIP
/
ITRPAPSM.ZIP
/
DAPLL.FOR
< prev
next >
Wrap
Text File
|
1985-11-29
|
5KB
|
125 lines
C
C ..................................................................
C
C SUBROUTINE DAPLL
C PURPOSE
C SET UP NORMAL EQUATIONS FOR A LINEAR LEAST SQUARES FIT
C TO A GIVEN DISCRETE FUNCTION
C
C USAGE
C CALL DAPLL(FFCT,N,IP,P,WORK,DATI,IER)
C SUBROUTINE FFCT REQUIRES AN EXTERNAL STATEMENT
C
C DESCRIPTION OF PARAMETERS
C FFCT - USER CODED SUBROUTINE WHICH MUST BE DECLARED
C EXTERNAL IN THE MAIN PROGRAM. IT IS CALLED
C CALL FFCT(I,N,IP,P,DATI,WGT,IER) AND RETURNS
C THE VALUES OF THE FUNDAMENTAL FUNCTIONS FOR
C THE I-TH ARGUMENT IN P(1) UP TO P(IP)
C FOLLOWED BY THE I-TH FUNCTION VALUE IN P(IP+1)
C N IS THE NUMBER OF ALL POINTS
C P,DATI,WGT MUST BE OF DOUBLE PRECISION.
C DATI IS A DUMMY PARAMETER WHICH IS USED AS ARRAY
C NAME. THE GIVEN DATA SET MAY BE ALLOCATED IN DATI
C WGT IS THE WEIGHT FACTOR FOR THE I-TH POINT
C IER IS USED AS RESULTANT ERROR PARAMETER IN FFCT
C N - NUMBER OF GIVEN POINTS
C IP - NUMBER OF FUNDAMENTAL FUNCTIONS USED FOR LEAST
C SQUARES FIT
C IP SHOULD NOT EXCEED N
C P - WORKING STORAGE OF DIMENSION IP+1, WHICH
C IS USED AS INTERFACE BETWEEN APLL AND THE USER
C CODED SUBROUTINE FFCT
C P MUST BE OF DOUBLE PRECISION.
C WORK - WORKING STORAGE OF DIMENSION (IP+1)*(IP+2)/2.
C ON RETURN WORK CONTAINS THE SYMMETRIC COEFFICIENT
C MATRIX OF THE NORMAL EQUATIONS IN COMPRESSED FORM,
C I.E. UPPER TRINGULAR PART ONLY STORED COLUMNWISE.
C THE FOLLOWING IP POSITIONS CONTAIN THE RIGHT
C HAND SIDE AND WORK((IP+1)*(IP+2)/2) CONTAINS
C THE WEIGHTED SQUARE SUM OF THE FUNCTION VALUES
C WORK MUST BE OF DOUBLE PRECISION.
C DATI - DUMMY ENTRY TO COMMUNICATE AN ARRAY NAME BETWEEN
C MAIN LINE AND SUBROUTINE FFCT.
C DATI MUST BE OF DOUBLE PRECISION.
C IER - RESULTING ERROR PARAMETER
C IER =-1 MEANS FORMAL ERRORS IN SPECIFIED DIMENSIONS
C IER = 0 MEANS NO ERRORS
C IER = 1 MEANS ERROR IN EXTERNAL SUBROUTINE FFCT
C
C REMARKS
C TO ALLOW FOR EASY COMMUNICATION OF INTEGER VALUES
C BETWEEN MAINLINE AND EXTERNAL SUBROUTINE FFCT, THE ERROR
C PARAMETER IER IS TREATED AS A VECTOR OF DIMENSION 1 WITHIN
C SUBROUTINE DAPLL. ADDITIONAL COMPONENTS OF IER MAY BE
C INTRODUCED BY THE USER FOR COMMUNICATION BACK AND FORTH.
C IN THIS CASE, HOWEVER, THE USER MUST SPECIFY IER AS A
C VECTOR IN HIS MAINLINE.
C EXECUTION OF SUBROUTINE DAPLL IS A PREPARATORY STEP FOR
C CALCULATION OF THE LINEAR LEAST SQUARES FIT.
C NORMALLY IT IS FOLLOWED BY EXECUTION OF SUBROUTINE DAPFS
C
C SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
C THE EXTERNAL SUBROUTINE FFCT MUST BE FURNISHED BY THE USER
C
C METHOD
C HANDLING OF THE GIVEN DATA SET (ARGUMENTS,FUNCTION VALUES
C AND WEIGHTS) IS COMPLETELY LEFT TO THE USER
C ESSENTIALLY HE HAS THREE CHOICES
C (1) THE I-TH VALUES OF ARGUMENT, FUNCTION VALUE AND WEIGHT
C ARE CALCULATED WITHIN SUBROUTINE FFCT FOR GIVEN I.
C (2) THE I-TH VALUES OF ARGUMENT, FUNCTION VALUE AND WEIGHT
C ARE DETERMINED BY TABLE LOOK UP. THE STORAGE LOCATIONS
C REQUIRED ARE ALLOCATED WITHIN THE DUMMY ARRAY DATI
C (POSSIBLY IN P TOO, IN EXCESS OF THE SPECIFIED IP + 1
C LOCATIONS).
C ANOTHER POSSIBILITY WOULD BE TO USE COMMON AS INTERFACE
C BETWEEN MAIN LINE AND SUBROUTINE FFCT AND TO ALLOCATE
C STORAGE FOR THE DATA SET IN COMMON.
C (3) THE I-TH VALUES OF ARGUMENT, FUNCTION VALUE AND WEIGHT
C ARE READ IN FROM AN EXTERNAL DEVICE. THIS MAY BE EASILY
C ACCOMPLISHED SINCE I IS USED STRICTLY INCREASING FROM
C ONE UP TO N WITHIN APLL
C
C ..................................................................
C
SUBROUTINE DAPLL(FFCT,N,IP,P,WORK,DATI,IER)
C
C
C DIMENSIONED DUMMY VARIABLES
DIMENSION P(1),WORK(1),DATI(1),IER(1)
DOUBLE PRECISION P,WORK,DATI,WGT,AUX
C
C CHECK FOR FORMAL ERRORS IN SPECIFIED DIMENSIONS
IF(N)10,10,1
1 IF(IP)10,10,2
2 IF(N-IP)10,3,3
C
C SET WORKING STORAGE AND RIGHT HAND SIDE TO ZERO
3 IPP1=IP+1
M=IPP1*(IP+2)/2
IER(1)=0
DO 4 I=1,M
4 WORK(I)=0.D0
C
C START GREAT LOOP OVER ALL GIVEN POINTS
DO 8 I=1,N
CALL FFCT(I,N,IP,P,DATI,WGT,IER)
IF(IER(1))9,5,9
5 J=0
DO 7 K=1,IPP1
AUX=P(K)*WGT
DO 6 L=1,K
J=J+1
6 WORK(J)=WORK(J)+P(L)*AUX
7 CONTINUE
8 CONTINUE
C
C NORMAL RETURN
9 RETURN
C
C ERROR RETURN IN CASE OF FORMAL ERRORS
10 IER(1)=-1
RETURN
END