Google
 

Trailing-Edge - PDP-10 Archives - decuslib10-01 - 43,50144/curfit.bas
There are 2 other files named curfit.bas in the archive. Click here to see a list.
01000 DATA 5E37,5E37
01001 DATA 0,0
01002 PRINT "DATE OF LAST REVISION ";
01003 PRINT "07/22/68."
01004 PRINT
01005 READ P
01010 DIM X(200),Y(200),U(200),V(200)
01015 IF P=5E37 THEN 9995
01016 FORK=1TO6
01017LETF(K)=1
01018NEXTK
01020 RESTORE
01021PRINT
01022PRINT"PLEASE SPECIFY THE NUMBER OF VALUES (N) GIVEN AS DATA"
01023PRINT"FOR THE TWO INPUT VARIABLES, AND THE OUTPUT CODE (D)."
01024PRINT"(D=1 IF OUTPUT IS TO BE IN ORDER OF INCREASING VALUES"
01025PRINT"OF THE INDEPENDENT VARIABLE, ELSE D=0).  N,D = ";
01027 INPUT N,D
01028 PRINT
01030FORI=1TON
01031READY(I)
01032NEXTI
01033FORI=1TON
01034READX(I)
01035NEXTI
01036READX,Y
01037IFX<5E37THEN9997
01038IFY<5E37THEN9997
01040PRINT
01045PRINT
01055 PRINT "       ";
01065 PRINT "L E A S T   S Q U A R E S   C U R V E   F I T"
01070 PRINT
01075 PRINT "CURVE TYPE","  INDEX  OF","   A","   B"
01080 PRINT " ","DETERMINATION"
01085 PRINT
01090 FOR I =1 TO 6
01095 FOR I1=1 TO 6
01100 LET S(I1)=0
01105 NEXT I1
01110 GOSUB 7000
01115 IF (I-5)*(I-6)=0 THEN 1195
01120 IF (I-2)*(I-3)=0 THEN 1155
01125 FOR J = 1 TO N
01130 LET V(J)=Y(J)
01135 GOSUB 5000
01140 NEXT J
01145 IF I=1 THEN 1245
01150 GO TO 1300
01155 FOR J = 1 TO N
01160 IF Y(J)<=0 THEN 1230
01165 LET V(J)=LOG(Y(J))
01170 GOSUB 5000
01175 NEXT J
01180 IF I=3 THEN 1270
01185 GO TO 1245
01195 FOR J = 1 TO N
01200 IF Y(J)=0 THEN 1230
01205 LET V(J)=1/Y(J)
01210 GOSUB 5000
01215 NEXT J
01220 IF I=6 THEN 1300
01225 GO TO 1245
01230 PRINT "DATA DOES NOT CONFORM TO THIS GENERAL CURVE TYPE"
01235 PRINT " ","     BEGINNING WITH POINT ";J
01236LETF(I)=0
01240 GO TO 1335
01245 FOR J = 1 TO  N
01250 LET U(J)=X(J)
01255 GOSUB 5015
01260 NEXT J
01265 GO TO 1325
01270 FOR J = 1 TO N
01275 IF X(J)<=0 THEN 1230
01280 LET U(J)=LOG(U(J))
01285 GOSUB 5015
01290 NEXT J
01295 GO TO 1325
01300 FOR J = 1 TO N
01305 IF X(J)=0 THEN 1230
01310 LET U(J)=1/X(J)
01315 GOSUB 5015
01320 NEXT J
01325 GOSUB 8000
01330 PRINT C(I),A(I),B(I)
01335 NEXT I
01336 IFD<>1THEN1340
01338GOSUB6000
01340 PRINT
01342PRINT
01343PRINT
01345 PRINT "ENTER THE NUMBER OF THE CURVE FOR WHICH DETAILS"
01350 PRINT "ARE DESIRED OR TYPE AN 'S' TO STOP";
01355 INPUTI
01356LETK=I
01357LETD1=D
01360 PRINT
01361IFF(I)=1THEN1365
01362 GOSUB7000
01363PRINT " COULD NOT BE FIT."
01364GOTO1340
01365 GOSUB 7018
01370 IF(I-1)*(I-5)*(I-6)<>0THEN1420
01375FORJ=1TON
01380LETY=A(I)+B(I)*X(J)
01385IFI=1THEN1405
01390LETY=1/Y
01395IFI=5THEN1405
01400	LETY=X(J)/(A(I)*X(J)+B(I))
01405GOSUB9000
01410NEXTJ
01414LETD=D1
01415GOTO1340
01420FORJ=1TON
01425IFI=2THEN1455
01430IFI=3THEN1445
01435LETY=A(4)+B(4)/X(J)
01440GOTO1460
01445LETY=A(3)*(X(J)^B(3))
01450GOTO1460
01455LETY=A(2)*EXP(B(2)*X(J))
01460GOSUB9000
01465NEXTJ
01466 LETD=D1
01470GOTO1340
05000 LET S(5)=S(5)+V(J)^2
05005 LET S(3)=S(3)+V(J)
05010 RETURN
05015 LET S(1)=S(1)+U(J)
05020 LET S(2)=S(2)+U(J)^2
05025 LET S(4)=S(4)+U(J)*V(J)
05030 RETURN
06000 FORI=1TON-1
06010LETM=I
06020FORJ=I+1TON
06030IFX(M)<=X(J)THEN6050
06040LETM=J
06050NEXTJ
06060IFM=ITHEN6100
06065LETP=X(M)
06066LETQ=Y(M)
06070LETX(M)=X(I)
06071LETY(M)=Y(I)
06080LETX(I)=P
06081 LETY(I)=Q
06100NEXTI
06110RETURN
07000 LET K=I
07001IFK=1THEN 7016
07002IFK=2THEN 7014
07003IFK=3THEN 7012
07004IFK=4THEN 7010
07005IFK=5THEN 7008
07006 PRINT"6. Y=X/(A*X+B) ";
07007	RETURN
07008 PRINT "5. Y=1/(A+B*X) ";
07009	RETURN
07010PRINT"4. Y=A+(B/X)",
07011	RETURN
07012PRINT"3. Y=A*(X^B)",
07013	RETURN
07014PRINT"2. Y=A*EXP(B*X)";
07015	RETURN
07016PRINT"1. Y=A+(B*X)",
07017	RETURN
07018PRINT"   ";
07019GOSUB 7001
07020PRINT" IS A";
07021IFK=1THEN 7026
07022IFK=2THEN 7028
07023IFK=3THEN 7030
07024PRINT" HYPERBOLIC";
07025GOTO 7031
07026PRINT" LINEAR";
07027GOTO 7031
07028PRINT"N EXPONENTIAL";
07029GOTO 7031
07030PRINT" POWER";
07031PRINT" FUNCTION.  THE RESULTS"
07032IFK=1THEN 7034
07033PRINT"      OF A LEAST-SQUARES FIT OF ITS LINEAR TRANSFORM"
07034IFD<>1THEN 7036
07035PRINT"      (SORTED IN ORDER OF ASCENDING VALUES OF X)"
07036PRINT"      ARE AS FOLLOWS:"
07037PRINT
07038PRINT"X-ACTUAL","Y-ACTUAL"," Y-CALC"," PCT DIFFER"
07039PRINT
07041	RETURN
07042PRINT
08000 LET B=(N*S(4)-S(1)*S(3))/(N*S(2)-(S(1)^2))
08005 LET A=(S(3)-B*S(1))/N
08010 LET S1= S(5)-(S(3)^2)/N
08015 LET S2= (B^2)*(S(2)-(S(1)^2)/N)
08020 LET C(I)=S2/S1
08025 IF (I-1)*(I-4)*(I-5) = 0 THEN 8060
08030 IF (I-2)*(I-3) = 0 THEN 8050
08035 LET A(6)=A
08040 LET B(6)=B
08045 RETURN
08050 LET A(I)=EXP(A)
08055 GO TO 8065
08060 LET A(I)=A
08065 LET B(I)=B
08070 RETURN
09000PRINTX(J),Y(J),Y,
09005 LET D=Y(J)-Y
09006 LET D=.1*SGN(D)*INT(1000*ABS(D/Y))
09010IFD<0THEN9020
09011 IF D>0 THEN 9015
09012 PRINT "      0"
09013	RETURN
09015PRINT"       ";
09020PRINTD
09025RETURN
09995 PRINT"LIST THE FILE 'STADES*' FOR INSTRUCTIONS."
09996 STOP
09997 PRINT"WRONG AMOUNT OF DATA SPECIFIED???   PLEASE CHECK."
09998 STOP
09999 END