|
|
DataMuseum.dkPresents historical artifacts from the history of: CP/M |
This is an automatic "excavation" of a thematic subset of
See our Wiki for more about CP/M Excavated with: AutoArchaeologist - Free & Open Source Software. |
top - metrics - download
Length: 14208 (0x3780)
Types: TextFile
Names: »KUFI.PAS«
└─⟦692ac107c⟧ Bits:30005923 PolyPascal-80 V3.10 arbejdsdiskette 1
└─⟦this⟧ »KUFI.PAS«
PROGRAM CURVEFIT;
LABEL 5,8,10,15,20,35,40,45,50,70,99,110,220,330,440;
CONST MAX=100;
SP=32;(*SPACE*)
FF=12;(*CLEAR SCREEN*)
ESC=27;
U=999999.0;
VAR (*Prefix "S" stands for SUMMA*)
SX,SX2,SX3,SX4,SLNX,SLNX2,SLNX3,SLNX4
,SY,SY2,SLNY,SLNY2,SXY,SX2Y,SXLNY,SX2LNY
,SYLNX,SYLNX2,SLNXLNY,SLNX2LNY
,A1,A2,A3,B1,B2,B3,C1,C2,C3,D1,D2,D3
,DENOM,A,B,C,SMIN,XINP,YCALC,SR2 ,XMIN
,YMIN,STDEV,R:REAL;
I,NP,NR,K,BEST,CALCXY,FIT,NLINES,OM,SAVEOM
,COPY,SWAP:INTEGER;
CH:STRINGÆ1Å;
X: ARRAYÆ1..MAXÅ OF REAL;
Y: ARRAYÆ1..MAXÅ OF REAL;
S: ARRAYÆ1..8Å OF REAL;
PROCEDURE TOPLINE(S:STRINGÆ48Å);
(*UNSCROLLED TEXT ON NASCOM SCREEN*)
CONST TOPLINEADR=$0BC9;
VAR P:INTEGER;
BEGIN
FOR P:=1 TO LENGTH(S) DO
MEMÆP+TOPLINEADRÅ:=ORD(MID(S,P,1));
FOR P:=P TO 48 DO
MEMÆP+TOPLINEADRÅ:=SP;
END;
PROCEDURE DELAY(D:INTEGER);
VAR I:INTEGER;
BEGIN FOR I:=1 TO D*1000 DO;
END;
PROCEDURE SPCONT;
LABEL 99;
BEGIN
IF OM>0 THEN GOTO 99;
SCREEN(12,14);WRITE('PRESS SPACE TO CONTINUE!');
REPEAT UNTIL KEYBOARD=SP;DELAY(8);
99:
END;
(*PRINTER ROUTINES*)
PROCEDURE PRT;
CODE $D7,$02,$75,$00,$E1,$DF,$71;
PROCEDURE CRTPRT;
CODE $D7,$03,$65,$75,$00,$E1,$DF,$71;
PROCEDURE CRT;
CODE $DF,$77;
PROCEDURE PRINTER;CODE
$D7,$50,$FE,$98,$20,$02,$3E,$95,$FE,$D2,
$20,$02,$3E,$81,$FE,$F8,$20,$02,$3E,$89,
$FE,$AB,$20,$02,$3E,$BA,$D7,$2C,$FE,$0D,
$C0,$3E,$0A,$D7,$25,$3A,$03,$08,$CB,$5F,
$C8,$3A,$00,$08,$3D,$32,$00,$08,$C0,$3A,
$02,$08,$32,$00,$08,$3A,$01,$08,$3C,$32,
$01,$08,$06,$0C,$3E,$0A,$D7,$04,$10,$FA,
$3E,$0D,$F5,$DB,$00,$17,$38,$FB,$F1,$DF,
$6F,$C9,$E1,$22,$78,$0C,$3E,$C3,$32,$77,
$0C;
PROCEDURE OMODE(X:INTEGER);
(*SETS OUTPUT MODE.*)
BEGIN
OM:=X;
CASE OM OF
0: CRT;
1: PRT;
2:CRTPRT;
END;
END;
PROCEDURE CS;(*CLEAR SCREEN*)
BEGIN SAVEOM:=OM;OMODE(0);WRITE(CHR(FF));OMODE(SAVEOM);
END;
PROCEDURE PP(PARAM:STRINGÆ4Å);
BEGIN
WRITE(CHR(27));(*ESC*)
FOR I:=1 TO LENGTH(PARAM) DO WRITE(MID(PARAM,I,1));
END;
PROCEDURE PRINTPARAM;
(*PROGRAMMING PRINTER*)
BEGIN
(*PAGE PARAMETERS*)
INIT MEMÆ$800Å TO $3C,$00,$3C;(*60 LINES/PAGE*)
MEMÆ$803Å:=$08;(*FORM FEEDS*)
(*MEMÆ$803Å:=$00 NO FORM FEED*)
PRINTER;OMODE(1);
WRITE(CHR(15));(*ELONGATED? 14=+,15=- *)
PP('F');
PP('N');
PP('A');
PP('"');
PP('Y');
PP('L010');
PP('&');
OMODE(0);
END;
PROCEDURE RSQRT(I:INTEGER);
BEGIN
WRITELN('SUM OF SQUARED RESIDUALS =',SR2);WRITELN;
SR2:=ABS(SR2); STDEV:=SQRT(SR2/(NP-3)); SÆIÅ:=STDEV;
END;
PROCEDURE EQUATIONS(TYPE:INTEGER);
(*SOLVES THREE LINEAR EQUATIONS IN THREE VARIABELS, A,B AND C*)
LABEL 99;
BEGIN
SAVEOM:=OM; CS;
DENOM:=A1*(B2*C3-B3*C2)+A2*(B3*C1-B1*C3)+A3*(B1*C2-B2*C1);
IF DENOM=0 THEN
BEGIN
WRITELN('DENOM =0!,SKIPPED');
REPEAT UNTIL KEYBOARD=SP;
GOTO 99;
END;
A:=(D1*(B3*C2-B2*C3)+D2*(B1*C3-B3*C1)+D3*(B2*C1-B1*C2))/DENOM;
B:=(D1*(A2*C3-A3*C2)+D2*(A3*C1-A1*C3)+D3*(A1*C2-A2*C1))/DENOM;
C:=(D1*(A3*B2-A2*B3)+D2*(A1*B3-A3*B1)+D3*(A2*B1-A1*B2))/DENOM;
CASE TYPE OF
1:BEGIN
WRITELN('1. QUADRATIC FIT:');
IF CALCXY=1 THEN OMODE(0);
SR2:=SY2-A*SX2Y-B*SXY-C*SY;
RSQRT(1+SWAP);
END;
2:BEGIN
WRITELN('2.EXPONENTIAL FIT:');
IF CALCXY=1 THEN OMODE(0);
SR2:=SLNY2-A*SX2LNY-B*SXLNY-C*SLNY;
RSQRT(2+SWAP);
END;
3:BEGIN
WRITELN('3.LOGARITHMIC FIT:');
IF CALCXY=1 THEN OMODE(0);
SR2:=SY2-A*SYLNX2-B*SYLNX-C*SY;
RSQRT(3+SWAP);
END;
4:BEGIN
WRITELN('4.POWER FIT:');
IF CALCXY=1 THEN OMODE(0);
SR2:=SLNY2-A*SLNX2LNY-B*SLNXLNY-C*SLNY;
RSQRT(4+SWAP);
END;
END;
IF STDEV=U THEN GOTO 99;
WRITELN('A =',A:16:8);
WRITELN('B =',B:16:8);
WRITELN('C =',C:16:8);
WRITELN;
WRITELN('STANDARD DEVIATION =',STDEV:10:6);
WRITELN;WRITELN;SPCONT;
99:OMODE(SAVEOM);
END;
PROCEDURE NL;
BEGIN
NLINES:=NLINES+1;
IF NLINES=13 THEN
BEGIN
NLINES:=0;WRITELN;
IF CH<>'P' THEN SPCONT;
END;
END;
PROCEDURE SORTPOINTS;
BEGIN
FOR I:=1 TO NP DO
BEGIN
XMIN:=XÆIÅ;YMIN:=YÆIÅ;
FOR K:=I+1 TO NP DO
IF XMIN>XÆKÅ THEN
BEGIN
XÆIÅ:=XÆKÅ;YÆIÅ:=YÆKÅ;
XÆKÅ:=XMIN;YÆKÅ:=YMIN;
XMIN:=XÆIÅ;YMIN:=YÆIÅ;
END;
END;
END;
PROCEDURE LIST;
LABEL 99;
BEGIN
CS;TOPLINE ('**** Curve fitting program ***');
CASE CH OF
'D','L':OMODE(0);
'P': OMODE(2);
'T': OMODE(4);
END;
SORTPOINTS;
IF NP=0 THEN GOTO 99;
WRITELN('OBSERVED X,Y-VALUES:');WRITELN;
WRITELN(' X Y');
NLINES:=3;
FOR I:=1 TO NP DO
BEGIN WRITELN(I:3,'. ',XÆIÅ:11:6,' ',YÆIÅ:11:6);NL;
END;
99:
END;
PROCEDURE EDIT;
LABEL 10,99;
BEGIN
10:
SAVEOM:=OM;OMODE(0);
WRITE(CHR(ESC)
,'DEL,LIST,FIT,PRINT,SWAPXY,QUIT? ');
READ(CH);WRITE(CH);DELAY(4);
IF CH='P' THEN COPY:=1;
IF COPY=1 THEN PRINTPARAM
ELSE COPY:=0;
OMODE(SAVEOM);
CASE CH OF
'D': BEGIN(*DELETE*)
WRITELN;WRITE('DEL.POINT NR. ');
NR:=0;READ(NR);
IF (NR<1) OR (NR>NP) THEN
BEGIN WRITE(CHR(ESC));GOTO 99;
END;
FOR K:=NR TO NP-1 DO
BEGIN XÆKÅ:=XÆK+1Å; YÆKÅ:=YÆK+1Å;
END;
NP:=NP-1;LIST;
END;
'S' :BEGIN
FOR I:=1 TO NP DO
BEGIN R:=XÆIÅ;XÆIÅ:=YÆIÅ;YÆIÅ:=R;
END;
IF SWAP=0 THEN SWAP:=4 ELSE SWAP:=0;
LIST;
END;
'L','P':LIST;
'F','Q':;
OTHERS:GOTO 10;
END;
99:
END;
(************ MAIN PROGRAM ************)
BEGIN
5:(*RESET VARIABLES MAY BE USED IN A SECOND RUN!):*)
SX:=0;SX2:=0;SX3:=0;SX4:=0;SLNX:=0;
SLNX2:=0;SLNX3:=0;SLNX4:=0;SY:=0;
SY2:=0;SLNY:=0;SLNY2:=0;SXY:=0;
SX2Y:=0;SXLNY:=0;SX2LNY:=0;SYLNX:=0;
SYLNX2:=0;SLNXLNY:=0;SLNX2LNY:=0;
FOR I:=1 TO 8 DO SÆIÅ:=U;CALCXY:=0;
(*INPUT CORRESPONDING VALUES OF OBS. X AND Y:*)
8:
CH:='L';OMODE(0);LIST;
10:
NP:=NP+1;WRITE(NP:3,'. ');
XÆNPÅ:=U;READ(XÆNPÅ);WRITE(' ');
IF XÆNPÅ=U THEN
BEGIN NP:=NP-1;EDIT;
CASE CH OF
'F':GOTO 20;
'Q':GOTO 99;
OTHERS:GOTO 10;
END;
END;
15:
SCREEN(20,0);YÆNPÅ:=U;READ(YÆNPÅ);
IF YÆNPÅ=U THEN GOTO 15;WRITELN;
IF NP<MAX THEN GOTO 10;
20:
(*CALCULATE SUMS*)
SAVEOM:=OM;OMODE(0);
WRITELN;WRITELN('CALCULATING!');
OMODE(SAVEOM);
FOR I:=1 TO NP DO
BEGIN
SX:=SX+XÆIÅ;
SX2:=SX2+SQR(XÆIÅ);
SX3:=SX3+SQR(XÆIÅ)*XÆIÅ;
SX4:=SX4+SQR(SQR(XÆIÅ));
SLNX:=SLNX+LN(XÆIÅ);
SLNX2:=SLNX2+SQR(LN(XÆIÅ));
SLNX3:=SLNX3+LN(XÆIÅ)*SQR(LN(XÆIÅ));
SLNX4:=SLNX4+SQR(SQR(LN(XÆIÅ)));
SY:=SY+YÆIÅ;
SY2:=SY2+SQR(YÆIÅ);
SLNY:=SLNY+LN(YÆIÅ);
SLNY2:=SLNY2+SQR(LN(YÆIÅ));
SXY:=SXY+XÆIÅ*YÆIÅ;
SX2Y:=SX2Y+SQR(XÆIÅ)*YÆIÅ;
SXLNY:=SXLNY+XÆIÅ*LN(YÆIÅ);
SX2LNY:=SX2LNY+SQR(XÆIÅ)*LN(YÆIÅ);
SYLNX:=SYLNX+YÆIÅ*LN(XÆIÅ);
SYLNX2:=SYLNX2+YÆIÅ*SQR(LN(XÆIÅ));
SLNXLNY:=SLNXLNY+LN(XÆIÅ)*LN(YÆIÅ);
SLNX2LNY:=SLNX2LNY+SQR(LN(XÆIÅ))*LN(YÆIÅ);
END;
(*OUTPUT CALCULATED SUMS:*)
SPCONT;CS;
IF COPY=1 THEN
BEGIN
SAVEOM:=OM;OMODE(1);
WRITELN;WRITELN;WRITELN;
OMODE(SAVEOM);
END;
WRITELN('CALCULATED SUMS:');
WRITELN(CHR(171),'X = ',SX:20:10);
WRITELN(CHR(171),'X^2 = ',SX2:20:10);
WRITELN(CHR(171),'X^3 = ',SX3:20:10);
WRITELN(CHR(171),'X^4 = ',SX4:20:10);
WRITELN(CHR(171),'LN(X) = ',SLNX:20:10);
WRITELN(CHR(171),'LN(X)^2 = ',SLNX2:20:10);
WRITELN(CHR(171),'LN(X)^3 = ',SLNX3:20:10);
WRITELN(CHR(171),'LN(X)^4 = ',SLNX4:20:10);
WRITELN(CHR(171),'Y = ',SY:20:10);
WRITELN(CHR(171),'Y^2 = ',SY2:20:10);
WRITELN(CHR(171),'LN(Y) ',SLNY:20:10);
WRITELN(CHR(171),'LN(Y)^2 = ',SLNY2:20:10);
SPCONT;CS;
WRITELN(CHR(171),' X*Y ',SXY:20:10);
WRITELN(CHR(171),' Y*X^2 ',SX2Y:20:10);
WRITELN(CHR(171),' X*LN(Y) ',SXLNY:20:10);
WRITELN(CHR(171),' X^2*LN(Y)',SX2LNY:20:10);
WRITELN(CHR(171),' Y*LN(X) ',SYLNX:20:10);
WRITELN(CHR(171),' Y*LN(X)^2',SYLNX2:20:10);
WRITELN(CHR(171),' LNX*LNY ',SLNXLNY:20:10);
WRITELN(CHR(171),' LNX2*LNY ',SLNX2LNY:20:10);
WRITELN;WRITELN('NUMBER OF POINTS= ',NP);
WRITELN;WRITELN;WRITELN;
WRITELN('FITS:');SPCONT;CS;
(*CALCULATE COEFFICIENTS IN LINAER EQUATIONS*)
110:(*1. Linear case :*)
A1:=SX4;A2:=SX3;A3:=SX2;
B1:=SX3;B2:=SX2;B3:=SX;
C1:=SX2;C2:=SX;C3:=NP;
D1:=-SX2Y;D2:=-SXY;D3:=-SY;
EQUATIONS(1);
IF CALCXY=1 THEN GOTO 35;
220:(*2. Exponential case:*)
A1:=SX4;A2:=SX3;A3:=SX2;
B1:=SX3;B2:=SX2;B3:=SX;
C1:=SX2;C2:=SX;C3:=NP;
D1:=-SX2LNY;D2:=-SXLNY;D3:=-SLNY;
EQUATIONS(2);
IF CALCXY=1 THEN GOTO 35;
330:(*3. Logarithmic case:*)
A1:=SLNX4;A2:=SLNX3;A3:=SLNX2;
B1:=SLNX3;B2:=SLNX2;B3:=SLNX;
C1:=SLNX2;C2:=SLNX;C3:=NP;
D1:=-SYLNX2;D2:=-SYLNX;D3:=-SY;
EQUATIONS(3);
IF CALCXY=1 THEN GOTO 35;
440:(*4. Power case:*)
A1:=SLNX4;A2:=SLNX3;A3:=SLNX2;
B1:=SLNX3;B2:=SLNX2;B3:=SLNX;
C1:=SLNX2;C2:=SLNX;C3:=NP;
D1:=-SLNX2LNY;D2:=-SLNXLNY;D3:=-SLNY;
EQUATIONS(4);
IF CALCXY=1 THEN GOTO 35;
(*FIND BEST FIT*)
SMIN:=SÆ1Å;BEST:=1;
FOR I:=2 TO 8 DO
IF SÆIÅ<SMIN THEN
BEGIN SMIN:=SÆIÅ;
BEST:=I;
END;
CS;
CASE BEST OF
1: WRITE('FIT 1 BEST! (QUADRATIC)');
2: WRITE('FIT 2 BEST (EXPONENTIAL)');
3: WRITE('FIT 3 BEST (LOGARITHMIC)');
4: WRITE('FIT 4 BEST (POWER)');
5: WRITE('FIT 5 BEST! (QUADRATIC)');
6: WRITE('FIT 6 BEST! (EXPONENTIAL)');
7: WRITE('FIT 7 BEST! (LOGARITHMIC)');
8: WRITE('FIT 8 BEST! (POWER)');
END;
WRITELN(' (ST.DEV = ',SÆBESTÅ:9:6,')');
CALCXY:=1;
50:
WRITELN;FIT:=0;SAVEOM:=OM;OMODE(0);
WRITE('FIT 1,2,3,4 OR CR? ');READLN(FIT);
OMODE(SAVEOM);WRITELN;
IF FIT=0 THEN
BEGIN;
CALCXY:=0;OMODE(0);(*PRINTER OFF*);
WRITE('REPEAT,QUIT?');
70:READ(CH);
CASE CH OF
'R': BEGIN WRITELN;GOTO 5;
END;
'Q':GOTO 99;
OTHERS:GOTO 70;
END;
END;
CASE FIT OF
1:GOTO 110;
2:GOTO 220;
3:GOTO 330;
4: GOTO 440;
END;
35:CS;WRITELN;
40:
CASE FIT OF
1: TOPLINE('Y =A*X^2 + B*C (1. QUADRATIC)');
2: TOPLINE('Y =EXP(A*X^2 +B*X +C) (2. EXPONENTIAL)');
3:TOPLINE ('Y = A*LNX^2 +B*LNX +C) (3. LOGARITHMIC');
4:TOPLINE ('Y =EXP(A*LN^2+B*LNX+C) (4. POWER');
END;
IF STDEV=0 THEN GOTO 50;
WRITE('CALCULATED X,Y-VALUES:');
IF SWAP=4 THEN WRITELN('( X,Y SWAPPED!)')
ELSE WRITELN;
45:
XINP:=U;WRITE('X: ');
SAVEOM:=OM;OMODE(0);READ(XINP);
IF XINP =U THEN GOTO 50;
IF COPY=1 THEN
BEGIN OMODE(1); WRITE(XINP:11:6);OMODE(SAVEOM);
END;
CASE FIT OF
1:(*QUADRATIC:*)
YCALC:=A*SQR(XINP)+B*XINP+C;
2:(*EXPONENTIAL:*)
YCALC:=EXP(A*SQR(XINP)+B*XINP+C);
3:(*LOGARITHMIC*)
YCALC:=A*LN(SQR(XINP))+B*LN(XINP)+C;
4:(*POWER:)
YCALC:=EXP(A*LN(SQR(XINP))+B*LN(XINP)+C);
END;
WRITELN(' Y: ',YCALC:11:6);
GOTO 45;
99:WRITELN;CS;OMODE(0);
END.
1: (*Linear:*)
YCALC:=A*SQR(XINP)+B*XINP+C;
2:(*Exponential:*)
YCALC:=EXP(A*SQR(XINP)+B*XINP+C);
3:(*Logarithmic:*)
YCALC:=A*LN(SQR(XINP))+B*LN(XINP)+C;
4: (*Power:*)
YCALC:=EXP(A*LN(SQR(XINP))+B*LN(XINP)+C);
END;
SCREEN(20,0);
WRITELN('Y: ',YCALC:14:10);
GOTO 40;
99:
END.
«eof»