Numerical Smoothing/Differentation


Language: QBasic.

Objective: In some occasions one must calculate the derivative dY/dX from empirical data, for example, the heating rate of a furnace. As empirical data can not always be described by a mathematical function, this program does the required differentation direct from this data, using a numerical/geometrical method. However, the derivative calculated from empirical data is extremely susceptible to small deviations due to measuring errors, for example. So, the values of the derivatives can show a significative dispersion. To correct this problem, this program includes a Pseudo-Fourier smoothing algorithm, that can be used to smooth previously the empirical data, reducing the magnitude of aleatory errors. So, this program can be used in three ways: to smooth, to diferentiate or to smooth-and-diferentiate empirical data!


***** Begin of Program Listing *****

10 REM ***
20 REM ***     NUMERICAL DIFFERENTATION/SMOOTHING OF EMPIRICAL DATA
30 REM ***
32 REM ***
33 REM ***  References:
34 REM ***
35 REM ***    - LE DUY, A. & ZAJIC, J.E. Biotechnology and Bioengineering,
36 REM ***      Vol. XV, 1973, 805-810.
37 REM ***
38 REM ***    - AUBANEL, E.E. & OLDHAM, K.B. Byte, Feb. 1985, 207-218.       
39 REM ***
40 REM ***
41 REM ***                       Antonio Augusto Gorni
50 REM ***
60 REM ***          - Applesoft Basic Version: June 23, 1988
70 REM ***          - QuickBasic Version:      October 14, 1994
80 REM ***
100 KEY OFF
110 B$ = CHR$(7): OPTION BASE 1: PI = 3.141592654#
120  DIM TT(350), CX(2, 350), CXA(350), DX(2, 350), ZB(350), ZC(350), ZO(350), ZN(350), ZM(350), TC(350), ZI(350), CR(350), VS(30), U(180), V(180), R(180), I(180)
130  ON ERROR GOTO 3890: AQ$ = "": E = 1: FLAGSMOOTH = 0: INITPLOT = 0: FLAGPRINTER = 0: FLAGTRANSF = 0: FLAGSAV = 0: FLAGORDEN = 0
140  VIEW PRINT 1 TO 24: CLS : BF$ = "NUMERICAL SMOOTHING/DIFFERENTATION OF EMPIRICAL DATA": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$
150  COLOR 7, 0: LOCATE 6
160  PRINT TAB(31); "<1> Data Input"
170  PRINT TAB(31); "<2> Data Correction"
180  PRINT TAB(31); "<3> Data Suppression"
190  PRINT TAB(31); "<4> Data Sorting"
200  PRINT TAB(31); "<5> Data Transforming"
210  PRINT TAB(31); "<6> Data Output"
220  PRINT TAB(31); "<7> Smoothing"
230  PRINT TAB(31); "<8> Diferentation"
240  PRINT TAB(31); "<9> Results Output"
250  PRINT TAB(31); "<10> Graphics Plotting"
260  PRINT TAB(31); "<11> End"
270  LOCATE 22
280  PRINT TAB(31); "Your Choice"; : COLOR 7, 0: INPUT ""; R: IF R < 1 OR R > 11 THEN PRINT B$: GOTO 270
290  ON R GOTO 300, 640, 720, 850, 950, 1170, 2370, 1560, 1850, 3910, 2320
300  CLS
310  BF$ = "DATA INPUT": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 4 TO 24
320 IF FLAGSAV = 0 THEN 350
330 PRINT B$: LOCATE 12: INPUT "Current Data Not Saved Yet! Do You Want  to Continue (Y/N)"; R$
340 IF R$ = "Y" OR R$ = "y" THEN CLS  ELSE IF R$ = "N" OR R$ = "n" THEN 140 ELSE 330
350  FLAGTRANSF = 0: FLAGORDEN = 0: E = 1: LOCATE 7
360  COLOR 7, 0: PRINT TAB(31); "<1> New Data"
370  LOCATE 12
380  PRINT TAB(31); "<2> More Data"
390  LOCATE 19
400  PRINT TAB(31); "Your Choice"; : INPUT ""; RE: IF RE < 1 OR RE > 2 THEN PRINT B$:  GOTO 390
410  IF RE = 1 THEN N = 0
420  VIEW PRINT 1 TO 24: CLS : BF$ = "INPUT OPTION": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
430  LOCATE 7
440  PRINT TAB(31); "<1> Keyboard"
450  LOCATE 12
460  PRINT TAB(31); "<2> Disk"
470  LOCATE 19
480  PRINT TAB(31); "Your Choice"; : INPUT ""; RE: IF RE < 1 OR RE > 2 THEN PRINT B$: GOTO 470
490  ON RE GOTO 500, 550
500  AQ$ = "": CLS : BF$ = "INPUT VIA KEYBOARD": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: LOCATE 10: VIEW PRINT 4 TO 24: PRINT "Now Enter the Data": PRINT "Enter END to Finish"
510 N = N + 1: PRINT : PRINT "Point # "; N: INPUT "X"; B1$
520  IF B1$ = "END" OR B2$ = "end" THEN N = N - 1: FLAGSAV = 1: GOTO 140
530  INPUT "Y"; B2$: IF B2$ = "END" OR B2$ = "end" THEN N = N - 1: FLAGSAV = 1: GOTO 140
540 TT(N) = VAL(B1$): CX(1, N) = VAL(B2$): GOTO 510
550 CLS : BF$ = "INPUT VIA DISK": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
560 LOCATE 10: VIEW PRINT 4 TO 24: LOCATE 7: INPUT "Input File Name (? to List Directory)"; AQ$
570 IF AQ$ = "?" THEN PRINT : INPUT "Path"; R$: PRINT : FILES R$ + "*.FOU" ELSE 590
580 LOCATE 22: INPUT "Press  to Continue...", R$: CLS : GOTO 560
590  AQ$ = AQ$ + ".FOU": LOCATE 20: COLOR 31, 0: PRINT "Reading "; AQ$: COLOR 7, 0
600 OPEN AQ$ FOR INPUT AS #1: I = 1
610 INPUT #1, TT(N + I): INPUT #1, CX(1, N + I)
620 IF NOT EOF(1) THEN I = I + 1: GOTO 610
630  CLOSE #1: N = N + I: COLOR 7, 0: FLAGSAV = 0: GOTO 140
640  CLS : BF$ = "DATA CORRECTION": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: PRINT : PRINT : VIEW PRINT 4 TO 24: LOCATE 6
650 FLAGTRANSF = 0: FLAGSAV = 1: FLAGORDEN = 0: AQ$ = ""
660 ID = 1: IF E = 1 THEN 700
670 CLS : LOCATE 7: PRINT TAB(31); "<1> Original Data": LOCATE 12: PRINT TAB(31); "<2> Smoothed Data"
680 LOCATE 19: PRINT TAB(31); "Your Choice"; : INPUT ""; ID:  IF ID < 1 OR ID > 2 THEN BEEP: GOTO 680
690 CLS : IF ID = 1 THEN E = 1
700  INPUT "Index of the Point (=0 to End)"; IC: IF IC = 0 THEN 140
710  PRINT : PRINT "X = "; TT(IC); TAB(20); : INPUT TT(IC): PRINT "Y = "; CX(ID, IC); TAB(20); : INPUT CX(ID, IC): PRINT : GOTO 700
720  CLS : BF$ = "DATA SUPPRESSION": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: LOCATE 10: VIEW PRINT 4 TO 24: CONT1 = 0
730 FLAGTRANSF = 0: FLAGSAV = 1: AQ$ = "": FLAGORDEN = 0
740 ID = 1: IF E = 1 THEN 780
750 CLS : LOCATE 7: PRINT TAB(31); "<1> Original Data": LOCATE 12: PRINT TAB(31); "<2> Smoothed Data"
760 LOCATE 19: PRINT TAB(31); "Your Choice"; : INPUT ""; ID:  IF ID < 1 OR ID > 2 THEN BEEP: GOTO 760
770 CLS : IF ID = 1 THEN E = 1
780  LOCATE 11: INPUT "Index of the Point (=0 to End)"; IC: IF IC = 0 THEN 820
790  LOCATE 13: PRINT TT(IC); TAB(20); CX(ID, IC): PRINT : PRINT : PRINT B$; : INPUT "Confirm (Y/N)! ", BF$: IF BF$ = "N" OR BF$ = "n" THEN CLS : GOTO 780
800  IF BF$ <> "Y" AND BF$ <> "y" THEN 790
810  CLS : CONT1 = CONT1 + 1: VS(CONT1) = IC: GOTO 780
820  LOCATE 16: COLOR 31, 0: PRINT "Making Suppressions...": COLOR 7, 0: FOR I = 1 TO CONT1 - 1: FOR J = I + 1 TO CONT: IF VS(J) > VS(I) THEN IA = VS(J): VS(J) = VS(I): VS(I) = IA
830  NEXT J, I
840  FOR J = 1 TO CONT1: FOR I = VS(J) TO N - 1: TT(I) = TT(I + 1): CX(ID, I) = CX(ID, I + 1): NEXT I: N = N - 1: NEXT J: CLS : GOTO 140
850 IF FLAGORDEN = 1 THEN 140
860 CLS : BF$ = "DATA SORTING": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
870 IF AQ$ <> "" THEN LOCATE 9: PRINT "Original File: "; AQ$ ELSE PRINT "Original Data Not Saved Yet!"
880 RO$ = "": RA$ = "": RT$ = "": IF FLAGORDEN = 0 THEN RO$ = "NO" ELSE RO$ = "YES"
890 IF FLAGTRANSF = 0 THEN RT$ = "NO" ELSE RT$ = "YES"
900 IF E = 1 THEN RA$ = "NO" ELSE RA$ = "YES"
910 LOCATE 11: PRINT "Data Status -> Sorting: "; RO$; " - Transformation: "; RT$; " - Smoothing: "; RA$
920 LOCATE 13: PRINT "Cogito, Ergo Sum!": LOCATE 15: PRINT "This Can Take Some Time..."
930 FOR I = 1 TO N - 1: FOR J = I + 1 TO N: IF TT(I) > TT(J) THEN SWAP TT(I), TT(J): SWAP CX(ID, I), CX(ID, J)
940 NEXT J, I: FLAGORDEN = 1: GOTO 140
950  CLS : BF$ = "DATA TRANSFORMATION": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
960 VIEW PRINT 4 TO 24
970 IF FLAGTRANSF = 0 THEN FLAGTRANSF = 1: GOTO 1000
980 PRINT B$: LOCATE 12: INPUT "Transformation Done! Continue"; R$
990 IF R$ = "Y" OR R$ = "y" THEN 1000 ELSE IF R$ = "N" OR R$ = "n" THEN 140 ELSE 980
1000  LOCATE 7: PRINT "Data transforming can be done through a routine programmed from the line 10000 on!": PRINT : PRINT
1010  PRINT : PRINT "If you will program it now, please remember that:"
1020 PRINT TAB(15); "-> TT(i)   is the array that contains the values of the variable ;"
1030 PRINT TAB(15); "-> CX(1,i) is the array that contains the original values of the variable ;"
1040 PRINT TAB(15); "-> CX(2,i) is the array that contains the smoothed values of the variable ;"
1050 PRINT TAB(15); "-> DX(1,i) is the array that contains the values of the original derivative ;"
1060 PRINT TAB(15); "-> DX(2,i) is the array that contains the values of the smoothed derivative ;"
1070 PRINT TAB(15); "-> N is the number of data points.": PRINT : PRINT : PRINT
1080 LOCATE 20: INPUT "The Function Is Programmed Now (Y/N)"; R$: IF R$ <> "N" AND R$ <> "n" AND R$ <> "Y" AND R$ <> "y" THEN PRINT B$: GOTO 1080
1090  IF R$ = "N" OR R$ = "n" THEN 140
1100 CLS : IF AQ$ <> "" THEN LOCATE 9: PRINT "Original File: "; AQ$ ELSE PRINT "Original Data Not Saved Yet!"
1110 RO$ = "": RA$ = "": RT$ = "": IF FLAGORDEN = 0 THEN RO$ = "NO" ELSE RO$ = "YES"
1120 IF FLAGTRANSF = 0 THEN RT$ = "NO" ELSE RT$ = "YES"
1130 IF E = 1 THEN RA$ = "NO" ELSE RA$ = "YES"
1140 LOCATE 11: PRINT "Data Status -> Sorting: "; RO$; " - Transformation: "; RT$; " - Smoothing: "; RA$
1150 LOCATE 13: PRINT "Cogito, Ergo Sum!": LOCATE 15: PRINT "This Can Take Some Time..."
1160 GOSUB 10000: GOTO 140
1170 CLS : COLOR 0, 7: LOCATE , 33: PRINT "DATA OUTPUT"
1180 LOCATE 7
1190 COLOR 7, 0: PRINT TAB(31); "<1> Screen"
1200 LOCATE 11
1210 PRINT TAB(31); "<2> Printer"
1220 LOCATE 15
1230 PRINT TAB(31); "<3> Disk"
1240 LOCATE 21
1250 PRINT TAB(31); "Your Choice"; : INPUT ""; RE: IF RE < 1 OR RE > 3 THEN PRINT B$: GOTO 1240
1260 ON RE GOTO 1270, 1370, 1480
1270 CLS : BF$ = "DATA OUTPUT": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: PRINT : PRINT : VIEW PRINT 4 TO 24
1280 IF AQ$ <> "" THEN LOCATE 11: PRINT "Original File: "; AQ$ ELSE PRINT "Original Data Not Saved Yet!"
1290 RO$ = "": RA$ = "": RT$ = "": IF FLAGORDEN = 0 THEN RO$ = "NO" ELSE RO$ = "YES"
1300 IF FLAGTRANSF = 0 THEN RT$ = "NO" ELSE RT$ = "YES"
1310 IF E = 1 THEN RA$ = "NO" ELSE RA$ = "YES"
1320 LOCATE 13: PRINT "Data Status -> Sorting: "; RO$; " - Transformation: "; RT$; " - Smoothing: "; RA$
1330 LOCATE 22: INPUT "Press  to Continue!", R$: CLS
1340 CL = 0: FOR I = 1 TO N: PRINT "Point #"; I: PRINT "X = "; TT(I); TAB(20); "Y = "; CX(1, I): PRINT : CL = CL + 3
1350 IF CL >= 17 AND I < N THEN PRINT : PRINT "Press  to Continue..."; : CL = 0: INPUT "", R$: CLS
1360 NEXT: PRINT : INPUT "Press  to Continue.", R$: CLS : GOTO 140
1370 CLS : BF$ = "DATA OUTPUT VIA PRINTER": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: LOCATE 8: VIEW PRINT 4 TO 24
1380 LOCATE 7: PRINT "Enter Data Identification Message."; : PRINT : LINE INPUT "", BF$: PRINT : PRINT : PRINT
1390 PRINT "Prepare Printer;": PRINT : PRINT "Mark Start of Report;": PRINT
1400 PRINT "Press  to Continue!": NP = 1
1410 INPUT "", R$
1420 LPRINT "Data for Numerical Smoothing/Differentation - Page no. "; NP: LPRINT : LPRINT BF$: LPRINT
1430 IF AQ$ <> "" THEN LPRINT "Data from the File "; AQ$,
1440 LPRINT DATE$, TIME$: LPRINT : LPRINT : LPRINT
1450 CL = 5: FOR I = 1 TO N: LPRINT "#"; I; TAB(19); TT(I); TAB(43); CX(1, I): CL = CL + 1
1460 IF CL >= 55 AND I < N THEN NP = NP + 1: LPRINT CHR$(12): LPRINT "Data for Numerical Smoothing/Differentation - Page no. "; NP: LPRINT BF$: LPRINT : CL = 0
1470 NEXT: LPRINT CHR$(12): CLS : GOTO 140
1480 CLS : BF$ = "DATA SAVING ON DISK": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: LOCATE 10: VIEW PRINT 4 TO 24
1490 PRINT "File Name (? to List Directory)"; : INPUT ""; AQ$
1500 IF AQ$ <> "?" THEN 1530
1510 INPUT "Path"; R$: PRINT : FILES R$ + "*.FOU"
1520 LOCATE 22: INPUT "Press  to Continue...", R$: CLS : GOTO 1490
1530 AQ$ = AQ$ + ".FOU": LOCATE 20: COLOR 31, 0: PRINT "Saving in "; AQ$: COLOR 7, 0
1540 OPEN AQ$ FOR OUTPUT AS #1: FOR I = 1 TO N: WRITE #1, TT(I): WRITE #1, CX(1, I): NEXT: CLOSE #1
1550 COLOR 7, 0: CLS : FLAGSAV = 0: GOTO 140
1560 VIEW PRINT 1 TO 24: CLS : BF$ = "NUMERICAL DIFFERENTATION": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7
1570 PRINT BF$: COLOR 7, 0
1580 VIEW PRINT 4 TO 24
1590 CLS : IF AQ$ <> "" THEN LOCATE 9: PRINT "Original File: "; AQ$ ELSE PRINT "Original Data Not Saved Yet!"
1600 RO$ = "": RA$ = "": RT$ = "": IF FLAGORDEN = 0 THEN RO$ = "NO" ELSE RO$ = "YES"
1610 IF FLAGTRANSF = 0 THEN RT$ = "NO" ELSE RT$ = "YES"
1620 IF E = 1 THEN RA$ = "NO" ELSE RA$ = "YES"
1630 LOCATE 11: PRINT "Data Status -> Sorting: "; RO$; " - Transformation: "; RT$; " - Smoothing: "; RA$
1640 LOCATE 13: PRINT "Cogito, Ergo Sum!"
1650 LOCATE 15: PRINT "This Can Take Some Time..."
1660 FOR I = 1 TO N: IF E = 1 THEN CXA(I) = CX(1, I): IND = 1 ELSE CXA(I) = CX(2, I): IND = 2
1670 NEXT I
1680  FOR I = 1 TO N - 2: IF CXA(I) = CXA(I + 1) THEN CXA(I + 1) = CXA(I) * .9999
1690  IF CXA(I) = CXA(I + 2) THEN CXA(I + 2) = CXA(I + 1) * .9999#
1700 IF TT(I) = TT(I + 1) THEN TT(I + 1) = TT(I) * .9999#
1710 IF TT(I) = TT(I + 2) THEN TT(I + 2) = TT(I + 1) * .9999#
1720  NEXT
1730 IF CXA(N) = CXA(N - 1) THEN CXA(N) = CXA(N) * .9999#
1740 IF TT(N) = TT(N - 1) THEN TT(N) = TT(N) * .9999#
1750  I = 1: DX(IND, I) = (CXA(I + 1) - CXA(I)) / (TT(I + 1) - TT(I))
1760 ZB(I + 1) = (CXA(I + 1) - CXA(I)) / (TT(I + 1) - TT(I)): ZC(I + 1) = (CXA(I + 2) - CXA(I + 1)) / (TT(I + 2) - TT(I + 1)): IF ABS(ZB(I + 1) - ZC(I + 1)) > .001 THEN 1810
1770 DX(IND, I + 1) = .5 * (ZB(I + 1) + ZC(I + 1))
1780 I = I + 1: IF I - N <= -2 THEN 1760
1790  IF I < N THEN 1840
1800 GOTO 140
1810 ZO(I + 1) = (TT(I + 1) - TT(I + 2)) / (CXA(I + 2) - CXA(I + 1)): ZN(I + 1) = .5 * (CXA(I + 1) + CXA(I + 2)) - ZO(I + 1) * .5 * (TT(I + 1) + TT(I + 2))
1820 ZM(I + 1) = (TT(I) - TT(I + 1)) / (CXA(I + 1) - CXA(I)): ZI(I + 1) = .5 * (CXA(I) + CXA(I + 1)) - ZM(I + 1) * .5 * (TT(I) + TT(I + 1))
1830 TC(I + 1) = (ZN(I + 1) - ZI(I + 1)) / (ZM(I + 1) - ZO(I + 1)): CR(I + 1) = ZM(I + 1) * TC(I + 1) + ZI(I + 1): DX(IND, I + 1) = (TC(I + 1) - TT(I + 1)) / (CXA(I + 1) - CR(I + 1)): GOTO 1780
1840 DX(IND, I + 1) = (CXA(I + 1) - CXA(I)) / (TT(I + 1) - TT(I)): GOTO 1780
1850  CLS : BF$ = "RESULTS OUTPUT": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
1860  LOCATE 7
1870  COLOR 7, 0: PRINT TAB(31); "<1> Screen"
1880  LOCATE 11
1890  PRINT TAB(31); "<2> Printer"
1900  LOCATE 15
1910  PRINT TAB(31); "<3> Disk"
1920  LOCATE 21
1930  PRINT TAB(31); "Your Choice"; : INPUT ""; RE: IF RE < 1 OR RE > 3 THEN PRINT B$: GOTO 1920
1940  ON RE GOTO 1950, 2080, 2230
1950  CLS : BF$ = "RESULTS OUTPUT": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: PRINT : PRINT : VIEW PRINT 4 TO 24
1960 IF AQ$ <> "" THEN LOCATE 10: PRINT "Original File: "; AQ$ ELSE PRINT "Original Data Not Saved Yet!"
1970 RO$ = "": RA$ = "": RT$ = "": IF FLAGORDEN = 0 THEN RO$ = "NO" ELSE RO$ = "YES"
1980 IF FLAGTRANSF = 0 THEN RT$ = "NO" ELSE RT$ = "YES"
1990 IF E = 1 THEN RA$ = "NO" ELSE RA$ = "YES"
2000 LOCATE 12: PRINT "Data Status -> Sorting: "; RO$; " - Transformation: "; RT$; " - Smoothing: "; RA$
2010 LOCATE 14: PRINT "Smoothing Degree: "; STR$(EPCT); "%": LOCATE 22: INPUT "Press  to Continue...", R$: CLS
2020 CL = 0: FOR I = 1 TO N: PRINT "Point #"; I: PRINT "X = "; TT(I); TAB(20);
2030 PRINT "Y = "; : IF E = 1 THEN PRINT CX(1, I) ELSE PRINT CX(2, I)
2040 PRINT "dX/dY = "; : IF E = 1 THEN PRINT DX(1, I) ELSE PRINT DX(2, I)
2050 PRINT : CL = CL + 4
2060 IF CL >= 16 AND I < N THEN PRINT : PRINT "Press  to Continue."; : CL = 0: INPUT "", R$: CLS
2070 NEXT: PRINT : INPUT "Press  to Continue.", R$: CLS : GOTO 140
2080 CLS : BF$ = "RESULTS OUTPUT VIA PRINTER": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: LOCATE 8: VIEW PRINT 4 TO 24
2090 LOCATE 7: PRINT "Enter Data Identification Message."; : PRINT : LINE INPUT "", BF$: PRINT : PRINT : PRINT
2100 PRINT B$; "Prepare Printer;": PRINT : PRINT "Mark Start of Report;": PRINT : PRINT "Press  to Continue.": NP = 1
2110 INPUT "", R$
2120 LPRINT "NUMERICAL SMOOTHING/DIFFERENTATION RESULTS - Page no. "; NP: LPRINT : LPRINT BF$: LPRINT
2130 IF AQ$ <> "" THEN LPRINT "Data from the File "; AQ$,
2140 LPRINT DATE$, TIME$: LPRINT
2150 LPRINT "Smoothing Degree: "; STR$(EPCT); "%": LPRINT : LPRINT : LPRINT : CL = 5
2160 FOR I = 1 TO N: LPRINT "#"; I; TAB(19); TT(I); TAB(37);
2170 IF E = 1 THEN LPRINT CX(1, I);  ELSE LPRINT CX(2, I);
2180 LPRINT TAB(55);
2190 IF E = 1 THEN LPRINT DX(1, I) ELSE LPRINT DX(2, I)
2200 CL = CL + 1
2210  IF CL >= 54 AND I < N THEN NP = NP + 1: LPRINT CHR$(12): LPRINT "NUMERICAL SMOOTHING/DIFFERENTATION RESULTS - Page no. "; NP: LPRINT BF$: LPRINT : CL = 0
2220  NEXT: LPRINT CHR$(12): VIEW PRINT 4 TO 24: CLS : GOTO 140
2230  CLS : BF$ = "RESULTS SAVING IN DISK": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: LOCATE 8: VIEW PRINT 4 TO 24
2240  CLS : LOCATE 8: PRINT "File Number (? to List Directory)"; : INPUT ""; AQ$
2250 IF AQ$ <> "?" THEN 2280
2260 PRINT : INPUT "Path"; R$: PRINT : R$ = R$ + "*.FOU": FILES R$
2270 PRINT : PRINT : INPUT "Press  to Continue...", R$: CLS : GOTO 2240
2280 AQ$ = AQ$ + ".FOU": LOCATE 22: COLOR 31, 0: PRINT "Saving in "; A1$: COLOR 7, 0
2290 OPEN AQ$ FOR OUTPUT AS #1: FOR I = 1 TO N: IF E = 1 THEN WRITE #1, TT(I), CX(1, I), DX(1, I) ELSE WRITE #1, TT(I), CX(2, I), DX(2, I)
2300 NEXT I: CLOSE #1
2310 GOTO 140
2320 VIEW PRINT 1 TO 24: IF FLAGSAV = 0 THEN 2360
2330 CLS : BF$ = "END OF PROGRAM RUN": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
2340 PRINT B$: LOCATE 12: INPUT "Current Data Not Saved! Do You Really Want to Quit"; R$
2350 IF R$ = "Y" OR R$ = "y" THEN 2360 ELSE IF R$ = "N" OR R$ = "n" THEN 140 ELSE GOTO 2340
2360 CLS : KEY ON: RESET: ON ERROR GOTO 0: END
2370 VIEW PRINT 1 TO 24: CLS : BF$ = "TREATMENT OF THE CURVE EXTREMES FOR SMOOTHING": PRINT TAB((80 - LEN(BF$)) / 2 + 1);
2380 COLOR 0, 7: PRINT BF$: COLOR 7, 0
2390 LOCATE 6
2400 COLOR 7, 0: PRINT TAB(31); "<1> No Correction of the Extremes"
2410 LOCATE 10
2420 PRINT TAB(31); "<2> Correction Using Straigth Line"
2430 LOCATE 14
2440 PRINT TAB(31); "<3> Mirror Effect"
2450 LOCATE 21
2460 PRINT TAB(31); "Your Choice"; : INPUT ""; RE: IF RE < 1 OR RE > 3 THEN BEEP: GOTO 2460
2470 FLAGEXTR = RE: ON RE GOTO 2680, 2590, 2480
2480 FOR I = 1 TO N: CXA(I) = CX(1, I): NEXT
2490 CLS : BF$ = "MIRROR POSITION ": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
2500 LOCATE 7, 27: PRINT "<1> Beginning from the Initial Point"
2510 LOCATE 12, 27: PRINT "<2> Beginning from the Final Point"
2520 LOCATE 21, 31: INPUT "Your Choice"; RE: IF RE < 1 OR RE > 3 THEN BEEP: GOTO 2520
2530 FLAGMIRROR = RE: ON RE GOTO 2540, 2570
2540 FOR I = N + 1 TO 2 * N: CXA(I) = CXA(I - N): NEXT
2550 FOR I = 1 TO N: CXA(I) = CX(1, N - I + 1): NEXT
2560 GOTO 2580
2570 FOR I = N + 1 TO 2 * N: CXA(I) = CX(1, 2 * N - I + 1): NEXT
2580 N = 2 * N: GOTO 2690
2590 CLS : BF$ = "NUMBER OF REFERENCE POINTS": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: FLAGRETA = 0
2600 VIEW PRINT 4 TO 24
2610 LOCATE 7, 31: PRINT "<1> Standard"
2620 LOCATE 12, 31: PRINT "<2> Defined by the User"
2630 LOCATE 21, 31: INPUT "Your Choice"; RE: IF RE < 1 OR RE > 3 THEN BEEP: GOTO 2520
2640 FLAGRETA = RE: ON RE GOTO 2680, 2650
2650 CLS : D1 = INT(N / 10): D2 = D1
2660 LOCATE 11: PRINT "No. of Reference Points, Lowest Extrem (Default ="; STR$(INT(N / 10)); ")"; : INPUT BF$: IF BF$ <> "" THEN D1 = VAL(BF$)
2670 LOCATE 13: PRINT "No. of Reference Points, Highest Extrem (Default ="; STR$(INT(N / 10)); ")"; : INPUT BF$: IF BF$ <> "" THEN D2 = VAL(BF$)
2680 FOR I = 1 TO N: CXA(I) = CX(1, I): NEXT
2690 VIEW PRINT 1 TO 24: CLS : BF$ = "NUMERICAL SMOOTHING": PRINT TAB((80 - LEN(BF$)) / 2 + 1);
2700 COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 4 TO 24
2710 Q = 0: N2 = INT((N + 1) / 2 + 1)
2720 IF AQ$ <> "" THEN LOCATE 8: PRINT "Original File: "; AQ$ ELSE PRINT "Original Data Not Saved Yet!"
2730 RO$ = "": RA$ = "": RT$ = "": IF FLAGORDEN = 0 THEN RO$ = "NO" ELSE RO$ = "YES"
2740 IF FLAGTRANSF = 0 THEN RT$ = "NO" ELSE RT$ = "YES"
2750 IF E = 1 THEN RA$ = "NO" ELSE RA$ = "YES"
2760 LOCATE 10: PRINT "Data Status -> Sorting: "; RO$; " - Transformation: "; RT$; " - Smoothing: "; RA$
2770 LOCATE 12:  INPUT "Smoothing Degree [%]"; EPCT
2780 E = INT((100 - EPCT) * (N + 1) / 200)
2790 IF E > INT((N + 1) / 2) THEN LOCATE 18: PRINT B$; "Excessive Degree! Decrease it!": LOCATE 22: INPUT "Press  to Continue...", R$
2800 IF E <> INT(E) OR E <= 1 THEN LOCATE 19: PRINT B$; "Incorrect Value! It must be Integer, Greater than One!": LOCATE 22: INPUT "Press  to Continue...", R$
2810 LOCATE 14: PRINT "Cogito, Ergo Sum!"
2820 LOCATE 16: PRINT "This Can Take Some Time..."
2830 IF E <= Q THEN 3530
2840 IF Q <> 0 THEN 3040
2850 IF FLAGEXTR = 1 OR FLAGEXTR = 3 THEN M = 0: B = 0: GOTO 2970
2860 IF FLAGRETA = 1 THEN D1 = INT(N / 10): D2 = D1
2870 S1 = 0: S2 = 0
2880 FOR J = 0 TO D1 - 1
2890 S1 = S1 + CXA(J + 1)
2900 NEXT J
2910 FOR J = 0 TO D2 - J
2920 S2 = S2 + CXA(N - J)
2930 NEXT J
2940 X1 = S1 / D1: X2 = S2 / D2
2950 M = (X2 - X1) / (N - (D1 + D2) / 2)
2960 B = (X1 + X2) / 2 - M * N / 2
2970 G = 0
2980 FOR J = 0 TO N - 1
2990 CXA(J + 1) = CXA(J + 1) - M * J - B
3000 G = G + CXA(J + 1)
3010 NEXT J
3020 R(1) = G / N
3030 Q = 1
3040 J2 = INT((N - 1) / 2)
3050 P1 = INT(LOG(2 * J2 - 1) / LOG(2))
3060  FOR K = Q TO E - 1
3070 J1 = J2
3080 S = PI * K * 2 / N
3090 C = COS(S): S = SIN(S)
3100  FOR J = 1 TO J1
3110 L = 2 * J - 1
3120 U(J + 1) = CXA(L + 1) * C + CXA(L + 2)
3130 V(J + 1) = CXA(L + 1) * S
3140  NEXT J
3150 S = 2 * S * C: C = 2 * C * C - 1
3160  FOR P = 1 TO P1
3170 U(J1 + 2) = 0: V(J1 + 2) = 0
3180 J1 = INT((J1 + 1) / 2)
3190  FOR J = 1 TO J1
3200 L = 2 * J - 1
3210 U = U(L + 1) * C - V(L + 1) * S + U(L + 2)
3220 V(J + 1) = U(L + 1) * S + V(L + 1) * C + V(L + 2)
3230 U(J + 1) = U
3240  NEXT J
3250 S = 2 * S * C: C = 2 * C * C - 1
3260  NEXT P
3270 R(K + 1) = (CXA(1) + (U(2) * C + V(2) * S)) / N
3280  NEXT K
3290  FOR K = Q TO E - 1
3300 J1 = J2
3310 S = 2 * PI * K / N
3320 C = COS(S): S = SIN(S)
3330  FOR J = 1 TO J1
3340 L = 2 * J - 1
3350 U(J + 1) = -(CXA(L + 1) * S)
3360 V(J + 1) = CXA(L + 1) * C + CXA(L + 2)
3370  NEXT J
3380 S = 2 * S * C: C = 2 * C * C - 1
3390  FOR P = 1 TO P1
3400 U(J1 + 2) = 0: V(J1 + 2) = 0
3410 J1 = INT((J1 + 1) / 2)
3420  FOR J = 1 TO J1
3430 L = 2 * J - 1
3440 U = U(L + 1) * C - V(L + 1) * S + U(L + 2)
3450 V(J + 1) = U(L + 1) * S + V(L + 1) * C + V(L + 2)
3460 U(J + 1) = U
3470  NEXT J
3480 S = 2 * S * C: C = 2 * C * C - 1
3490  NEXT P
3500 I(K + 1) = -((U(2) * C + V(2) * S) / N)
3510  NEXT K
3520  IF E > Q THEN Q = E
3530 F1 = 0: F2 = 0
3540  FOR K = 1 TO E - 1
3550 T = R(K + 1)
3560 F1 = F1 + T
3570 F2 = F2 + K * K * T
3580  NEXT K
3590 CX(2, 1) = R(1) + 2 * (F1 - F2 * (1 / E / E))
3600 CX(2, 1) = CX(2, 1) + B
3610 P1 = INT(LOG(2 * E - 3) / LOG(2))
3620  FOR J = 1 TO N - 1
3630 T2 = E * E
3640  FOR K = 1 TO E - 1
3650 F = 1 - K * K / T2
3660 U(K + 1) = R(K + 1) * F: V(K + 1) = -(I(K + 1) * F)
3670  NEXT K
3680 K1 = E - 1
3690 S = 2 * PI * J / N
3700 C = COS(S): S = SIN(S)
3710  FOR P = 1 TO P1
3720 U(K1 + 2) = 0: V(K1 + 2) = 0
3730 K1 = INT((K1 + 1) / 2)
3740  FOR K = 1 TO K1
3750 L = 2 * K - 1
3760 U = U(L + 1) * C - V(L + 1) * S + U(L + 2)
3770 V(K + 1) = U(L + 1) * S + V(L + 1) * C + V(L + 2)
3780 U(K + 1) = U
3790  NEXT K
3800 S = 2 * S * C: C = 2 * C * C - 1
3810  NEXT P
3820 CX(2, J + 1) = R(1) + 2 * (U(2) * C + V(2) * S)
3830 CX(2, J + 1) = CX(2, J + 1) + M * J + B
3840  NEXT J
3850 IF FLAGEXTR <> 3 THEN 140
3860 N = N / 2
3870 IF FLAGMIRROR = 1 THEN FOR I = 1 TO N: CX(2, I) = CX(2, I + N): NEXT I
3880 GOTO 140
3890  SCREEN 0: CLS : BEEP: LOCATE 12: COLOR 31, 0: PRINT "Error "; ERR; " at  line "; ERL; "!"
3900 COLOR 7, 0: LOCATE 21: INPUT "Press  to Continue!", BF$: RESUME 140
3910 IF INITPLOT = 0 THEN GRAFMODE$ = "POINT": INITPLOT = 1: DIM GRAPHICS(4141), LETRA$(255), X(350), Y(350): GOSUB 6200
3920 SM = 2: GOTO 4040
3930 X0 = X(1): X1 = X(1): Y0 = Y(1): Y1 = Y(1)
3940 FOR I = 2 TO N
3950 IF X0 > X(I) THEN X0 = X(I)
3960 IF X1 < X(I) THEN X1 = X(I)
3970 IF Y0 > Y(I) THEN Y0 = Y(I)
3980 IF Y1 < Y(I) THEN Y1 = Y(I)
3990 IF FLAGSMOOTH = 1 AND REVS = 1 THEN IF Y0 > CX(1, I) THEN Y0 = CX(1, I)
4000 IF FLAGSMOOTH = 1 AND REVS = 1 THEN IF Y1 < CX(1, I) THEN Y1 = CX(1, I)
4010 IF FLAGSMOOTH = 1 AND REVS = 2 THEN IF Y0 > DX(1, I) THEN Y0 = DX(1, I)
4020 IF FLAGSMOOTH = 1 AND REVS = 2 THEN IF Y1 < DX(1, I) THEN Y1 = DX(1, I)
4030 NEXT I: RETURN
4040  SCREEN 0: VIEW PRINT 1 TO 24: CLS : BF$ = "CURVE PLOTTING": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 2 TO 24
4050  LOCATE 6
4060  PRINT TAB(31); "<1> Define Graphic"
4070  PRINT : PRINT TAB(31); "<2> Define Axis"
4080  PRINT : PRINT TAB(31); "<3> Plot Graphic"
4090  PRINT
4100  PRINT TAB(31); "<4> Checker Graphic"
4110  PRINT
4120  PRINT TAB(31); "<5> Show Graphic": PRINT
4130  PRINT TAB(31); "<6> Print Graphic"
4140  PRINT
4150  PRINT TAB(31); "<7> Main Menu"
4160  LOCATE 23
4170  PRINT TAB(31); "Your Choice? "; : COLOR 7, 0: INPUT "", R: IF R < 1 OR R > 7 THEN BEEP: GOTO 4160
4180  ON R GOTO 4190, 4480, 4720, 6000, 6050, 6090, 140
4190 FLAGSMOOTH = 0: REVS = 0
4200  VIEW PRINT 1 TO 24: CLS : BF$ = "GRAPHIC DEFINITION": PRINT TAB((80 - LEN(BF$)) / 2 + 1);
4210  COLOR 0, 7: PRINT BF$: COLOR 7, 0
4220  LOCATE 7: PRINT TAB(31); "<1> X Axis"
4230  LOCATE 11
4240  PRINT TAB(31); "<2> Y Axis"
4250  LOCATE 15: PRINT TAB(31); "<3> Return to Graphics Menu": LOCATE 21
4260  PRINT TAB(31); "Your Choice"; : INPUT ""; REA: IF REA < 1 OR REA > 3 THEN PRINT B$: GOTO 4200
4270 IF REA = 1 THEN BF$ = "X-AXIS" ELSE IF REA = 2 THEN BF$ = "Y-AXIS" ELSE 4040
4280  CLS : BF$ = "CONTENT OF THE " + BF$: PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0
4290  VIEW PRINT 4 TO 24: LOCATE 7: PRINT TAB(31); "<1> Data, Original X-Axis "
4300  LOCATE 11
4310  PRINT TAB(31); "<2> Data, Original Y-Axis"
4320  LOCATE 15
4330  PRINT TAB(31); "<3> Values of dY/dX"
4340  LOCATE 21
4350  PRINT TAB(31); "Your Choice"; : INPUT ""; RE: IF RE < 1 OR RE > 3 THEN PRINT B$: GOTO 4280
4360  FOR I = 1 TO N
4370  IF RE = 1 THEN IF REA = 1 THEN X(I) = TT(I) ELSE Y(I) = TT(I)
4380  IF RE = 2 AND E = 1 THEN IF REA = 1 THEN X(I) = CX(1, I) ELSE Y(I) = CX(1, I)
4390  IF RE = 2 AND E <> 1 THEN REVS = 1: IF REA = 1 THEN X(I) = CX(2, I) ELSE Y(I) = CX(2, I)
4400  IF RE = 3 AND E = 1 THEN IF REA = 1 THEN X(I) = DX(1, I) ELSE Y(I) = DX(1, I)
4410  IF RE = 3 AND E <> 1 THEN REVS = 2: IF REA = 1 THEN X(I) = DX(2, I) ELSE Y(I) = DX(2, I)
4420  NEXT I
4430 IF RE = 1 OR REVS = 0 THEN 4470
4440 CLS
4450 LOCATE 12: INPUT "Do Your Want to Plot Simultaneously the Original and Smoothed Data"; R$
4460 IF R$ = "Y" OR R$ = "y" THEN FLAGSMOOTH = 1 ELSE IF R$ <> "N" AND R$ <> "n" THEN PRINT B$: GOTO 4450
4470 GOSUB 3930: VIEW PRINT 1 TO 24: GOTO 4200
4480 VIEW PRINT 1 TO 24: CLS : BF$ = "DEFINITION OF THE AXIS AND GRAPHIC PARAMETERS": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 4 TO 24
4490 LOCATE 7
4500 PRINT "Current  Minimum X: "; X0; : INPUT " - New: ", S$: IF S$ <> "" THEN X0 = VAL(S$)
4510 LOCATE 9
4520 PRINT "Current Maximum X: "; X1; : INPUT " - New: ", S$: IF S$ <> "" THEN X1 = VAL(S$)
4530 LOCATE 11
4540 PRINT "Current Minimum Y: "; Y0; : INPUT " - New: ", S$: IF S$ <> "" THEN Y0 = VAL(S$)
4550 LOCATE 13
4560 PRINT "Current Maximum Y: "; Y1; : INPUT " - New: ", S$: IF S$ <> "" THEN Y1 = VAL(S$)
4570 LOCATE 16
4580 PRINT "Current Number of Ticks in the X-Axis: "; NX: INPUT "New: ", S$: IF S$ <> "" THEN NX = VAL(S$)
4590 LOCATE 19
4600 PRINT "Current Number of Ticks in the Y-Axis: "; NY: INPUT "New: ", S$: IF S$ <> "" THEN NY = VAL(S$)
4610 CLS : LOCATE 7
4620 PRINT "Current X-Axis Label: "; XA$: LINE INPUT "New: ", S$: IF S$ <> "" THEN XA$ = S$
4630 LOCATE 10
4640 PRINT "Current Y-Axis Label: "; YA1$: LINE INPUT "New: ", S$: IF S$ <> "" THEN YA1$ = S$
4650 IF E <> 1 THEN YA$ = YA1$ + " - Smooth =" + STR$(EPCT) + "%" ELSE YA$ = YA1$
4660 LOCATE 13
4670 PRINT "Current Graphic Output: "; GRAFMODE$; : INPUT " - New (L/P/M): ", S$
4680 IF S$ <> "" THEN IF S$ = "P" OR S$ = "p" THEN GRAFMODE$ = "POINT" ELSE IF S$ = "L" OR S$ = "l" THEN GRAFMODE$ = "LINE" ELSE IF S$ = "P" THEN GRAFMODE$ = "POINT" ELSE IF S$ = "M" OR S$ = "m" THEN GRAFMODE$ = "MIX" ELSE 4660
4690  IF GRAFMODE$ = "LINE" THEN 4040 ELSE LOCATE 15
4700  PRINT "Current Point Type: "; SM; : INPUT " - Novo: ", S$: IF S$ <> "" THEN SM = VAL(S$)
4710  GOTO 4040
4720  VIEW PRINT 1 TO 24: CLS : BF$ = "GRAPHIC PLOTTING": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 4 TO 24
4730  LOCATE 11: PRINT "Press  to Plot Graphic.": PRINT : PRINT "After the BEEP the Graphic is Ready;": PRINT : PRINT "Then Press  to Continue.": LOCATE 23: INPUT "Now  to Plot! ", R$
4740  GOSUB 5300: GOSUB 5000
4750 IF GRAFMODE$ = "LINE" THEN 4830
4760  FOR I = 1 TO N
4770 X = X(I): Y = Y(I)
4780  GOSUB 5110
4790 IF PY < 12 OR PY > 177 THEN GOTO 4820
4800 IF PX < 51 OR PX > 639 THEN GOTO 4820
4810 IF SM = 1 THEN CIRCLE (PX, PY), 2.25 ELSE PSET (PX, PY), 0: DRAW LETRA$(SM + 26)
4820  NEXT: IF GRAFMODE$ <> "MIX" THEN 4910
4830 PTOINIC = 0: FOR I = 1 TO N
4840 X = X(I): Y = Y(I)
4850 GOSUB 5110
4860 IF PY < 12 OR PY > 177 THEN PTOINIC = 0: GOTO 4900
4870 IF PX < 51 OR PX > 639 THEN PTOINIC = 0: GOTO 4900
4880 IF PTOINIC = 0 THEN PTOINIC = 1: PX1 = PX: PY1 = PY: GOTO 4900
4890 LINE (PX1, PY1)-(PX, PY): PX1 = PX: PY1 = PY
4900 NEXT
4910 IF FLAGSMOOTH = 0 THEN 5910
4920  FOR I = 1 TO N
4930 X = X(I): IF REVS = 1 THEN Y = CX(1, I) ELSE Y = DX(1, I)
4940  GOSUB 5110
4950 IF PY < 12 OR PY > 177 THEN GOTO 4980
4960  IF PX < 51 OR PX > 639 THEN GOTO 4980
4970 CIRCLE (PX, PY), 2.25
4980  NEXT
4990  GOTO 5910
5000 IF EX = 0 THEN 5030 ELSE IF EX > 0 THEN OFS = 1 ELSE OFS = 0
5010 EX$ = "(x10^" + RIGHT$(STR$(3 * EX), LEN(STR$(3 * EX)) - OFS) + ")"
5020 X$ = X$ + " " + EX$
5030 IF EY = 0 THEN 5060 ELSE IF EY > 0 THEN OFS = 1 ELSE OFS = 0
5040 EY$ = "(x10^" + RIGHT$(STR$(3 * EY), LEN(STR$(3 * EY)) - OFS) + ")"
5050 Y$ = Y$ + " " + EY$
5060 X = 51 + INT(529 - 10 * LEN(X$)) / 2: Y = 198
5070 P$ = X$: GOSUB 5230
5080 X = 51: Y = 8
5090 P$ = Y$: GOSUB 5230
5100  RETURN
5110 PX = 51 + 588 * (X - X0) / (X1 - X0)
5120 PY = 12 + 165 * (Y1 - Y) / (Y1 - Y0)
5130  RETURN
5140 POV = 0
5150 SI = 1: IF VAR < 0 THEN SI = -1
5160 VAR = ABS(VAR): IF VAR < 10 ^ -(EXPO + 3) THEN VAR = 0: GOTO 5190
5170 IF VAR > 999 THEN VAR = VAR / 1000: POV = POV + 1: GOTO 5170
5180 IF VAR < 1 THEN VAR = VAR * 1000: POV = POV - 1: GOTO 5170
5190 VAR$ = STR$(SI * VAR)
5200 VAR$ = LEFT$(VAR$, FIXADOR)
5210 VAR = VAL(VAR$)
5220 RETURN
5230  FOR I = 1 TO LEN(P$)
5240 PA$ = MID$(P$, I, 1)
5250  IF PA$ = " " THEN 5280
5260 AP = ASC(PA$)
5270  PSET (X + (I - 1) * 10, Y), 0: DRAW LETRA$(AP)
5280  NEXT
5290  RETURN
5300 SCREEN 2: X$ = XA$: Y$ = YA$
5310 LINE (51, 12)-(639, 12): LINE (639, 12)-(639, 177): LINE (639, 177)-(51, 177): LINE (51, 177)-(51, 12)
5320 SX = (X1 - X0) / NX: SY = (Y1 - Y0) / NY
5330 FOR I = X0 TO X1 STEP SX
5340 X = I
5350 GOSUB 5110
5360 LINE (PX, 12)-(PX, 19): LINE (PX, 177)-(PX, 170)
5370 NEXT
5380 FOR I = Y0 TO Y1 STEP SY
5390 Y = I
5400 GOSUB 5110
5410 LINE (51, PY)-(67, PY): LINE (639, PY)-(622, PY)
5420 NEXT
5430 FIXADOR = 5
5440 EXPO = EX: VAR = X0: GOSUB 5140: PIV = POV
5450 VAR = X1: GOSUB 5140: IF ABS(PIV) < ABS(POV) THEN PIV = POV
5460 FOR K = X0 TO X1 STEP SX
5470 X = K: GOSUB 5110
5480 IF PX >= 603 THEN 5570
5490 VAR = K
5500 GOSUB 5140
5510 IF POV <> PIV THEN VAR = VAR * 10 ^ (3 * (POV - PIV))
5520 IF ABS(VAR) < .001 THEN VAR = 0
5530 IF PIV <> 0 THEN VAR = INT(VAR + .5)
5540 P$ = LEFT$(STR$(VAR), 5): Y = 187
5550 X = PX - 6 * LEN(P$)
5560 GOSUB 5230
5570 NEXT
5580 EX = PIV
5590 VAR = Y0: EXPO = EY: GOSUB 5140: PIV = POV
5600 VAR = Y1: GOSUB 5140: IF ABS(PIV) > ABS(POV) THEN PIV = POV
5610 FOR K = Y0 TO Y1 STEP SY
5620 Y = K: GOSUB 5110
5630 IF PY < 20 THEN 5720
5640 Y = PY + 2
5650 VAR = K
5660 GOSUB 5140
5670 IF POV <> PIV THEN VAR = VAR * 10 ^ (3 * (POV - PIV))
5680 IF ABS(VAR) < .001 THEN VAR = 0
5690 IF PIV <> 0 THEN VAR = INT(VAR + .5)
5700 P$ = LEFT$(STR$(VAR), 5): X = 51 - 10 * LEN(P$)
5710 GOSUB 5230
5720 NEXT
5730 EY = PIV
5740  IF X0 * X1 > 0 THEN 5770
5750 X = 0: GOSUB 5110
5760  LINE (PX, 12)-(PX, 177)
5770  IF Y0 * Y1 > 0 THEN 5800
5780 Y = 0: GOSUB 5110
5790  LINE (51, PY)-(639, PY)
5800  RETURN
5810 PU = 587 / NX: PA = 165 / NY
5820  FOR X = 51 TO 630 STEP PU
5830  FOR Y = 12 TO 177 STEP PA / 5
5840  PSET (X, Y)
5850  NEXT: NEXT
5860  FOR Y = 12 TO 177 STEP PA
5870  FOR X = 51 TO 639 STEP PU / 5
5880  PSET (X, Y)
5890  NEXT: NEXT
5900  RETURN
5910 GOTO 5980
5920  PIV = X0: FOR X = PIV TO X1 STEP (X1 - X0) / 588
5930  GOSUB 5110
5940  IF PY < 12 OR PY > 177 THEN 5970
5950  IF PX < 51 OR PX > 639 THEN 5970
5960  PSET (PX, PY), 1
5970  NEXT X
5980 GET (0, 0)-(639, 199), GRAPHICS
5990 A$ = INKEY$: IF A$ <> "" THEN 4040 ELSE 5990
6000  VIEW PRINT 1 TO 24: CLS : BF$ = "CHECKER GRAPHIC": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 4 TO 24
6010  LOCATE 11: PRINT "Now the Graphic Will Be Checkered.": PRINT : PRINT "After That, Press  to Continue.": LOCATE 23: INPUT "Now  to Checker!", R$
6020  SCREEN 2: PUT (0, 0), GRAPHICS, PSET: GOSUB 5810
6030 GET (0, 0)-(639, 199), GRAPHICS
6040 A$ = INKEY$: IF A$ <> "" THEN 4040 ELSE 6040
6050  VIEW PRINT 1 TO 24: CLS : BF$ = "SHOW GRAPHIC": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 4 TO 24
6060  LOCATE 11: PRINT "Press  to Show Graphic;": PRINT : PRINT "After That, Press  to Continue.": LOCATE 23: INPUT "Now,  to Show! ", R$
6070  SCREEN 2: PUT (0, 0), GRAPHICS, PSET
6080 A$ = INKEY$: IF A$ <> "" THEN 4040 ELSE 6040
6090  VIEW PRINT 1 TO 24: CLS : BF$ = "GRAPHIC PRINTING": PRINT TAB((80 - LEN(BF$)) / 2 + 1); : COLOR 0, 7: PRINT BF$: COLOR 7, 0: VIEW PRINT 2 TO 24
6100 LOCATE 11: PRINT "Prepare Printer; ": PRINT : INPUT "Then Press  to Continue! ", R$
6110  SCREEN 2: PUT (0, 0), GRAPHICS, PSET
6120 IF FLAGPRINTER = 0 THEN FLAGPRINTER = 1: RESET: OPEN "LPT1:" FOR RANDOM AS #2: WIDTH #2, 255
6130 PRINT #2, CHR$(24); : PRINT #2, CHR$(27); "A"; CHR$(8); : DEF SEG = &HB800
6140 FOR A = 0 TO 79: PRINT #2, CHR$(27); "K"; CHR$(144); CHR$(1); : B = A + &H1EF0
6150 FOR C = 1 TO 100: D = PEEK(B): EPRT = PEEK(B + &H2000): PRINT #2, CHR$(EPRT); CHR$(EPRT); CHR$(D); CHR$(D); : B = B - 80: NEXT
6160 PRINT #2, CHR$(13); CHR$(10); : NEXT
6170 PRINT #2, CHR$(13); CHR$(24); CHR$(27); CHR$(50);
6180 LPRINT CHR$(12)
6190 GOTO 4040
6200 NX = 5: NY = 5
6210 FOR I = 28 TO 125: READ LETRA$(I): NEXT
6220 FOR I = 1 TO 255: LETRA$(I) = "P1,0" + LETRA$(I): NEXT
6230 RETURN
6240 DATA "BL2D1R4U2L4D1"
6250 DATA "L4R8L4U2D4"
6260 DATA "BD2L2R4H3G3R1"
6270 DATA "BU2L2R4G3H3R1"
6280 DATA ""
6290 DATA "BR3U1BU2U3"
6300 DATA "BR1BU5U1BR3D1"
6310 DATA "BR2U6BR3D6U4R2L7BD2R7"
6320 DATA "BR4U6D1R2L4G1F1R4F1G1L4"
6330 DATA "BR1BU6D1R1U1BD6BR4U1L1D1BL4E6"
6340 DATA "BR6H5E1R2F1G4F1R2E2"
6350 DATA "BR3BU6D2"
6360 DATA "BR4H2U2E2"
6370 DATA "BR1E2U2H2"
6380 DATA "BU3R6L3U2D4U2H2F4H2E2G4"
6390 DATA "BU3R6L3U2D4"
6400 DATA "BR2E1U1L1"
6410 DATA "BU3R5"
6420 DATA "BR2U1R1D1"
6430 DATA "E6"
6440 DATA "R3E1U4H1L3G1D4E4"
6450 DATA "BR1R4L2U6G1"
6460 DATA "R5L5E5H1L3G1"
6470 DATA "BU1F1R3E1U1H2E2L5"
6480 DATA "BR4U6G4R5"
6490 DATA "BU1F1R3E1U2H1L4U2R5"
6500 DATA "BR1R3E1U1BD1G1L3BL1BU1U1BE1R3BL3BL1BE1E2"
6510 DATA "BR2U2E4L6"
6520 DATA "BR1R3E1U1H1L3BG1D1U1E1H1U1E1R3F1D1"
6530 DATA "BR1E4BU1H1L3G1D1F1R2BR1"
6540 DATA "BR2BU1U1R1D1L1BU3U1R1D1"
6550 DATA "BR2E1U1L1BU2U1R1D1"
6560 DATA "BU3F3H3E3"
6570 DATA "BU2R5BU2L5"
6580 DATA "BR5BU3H3F3G3"
6590 DATA "BR3U1BU2E2H1L3G1"
6600 DATA "BU1U1E1R1D3L1R3E1U4H1L3G1"
6610 DATA "U5BU1BR1R3BR1BD1D5BL5BU3R5"
6620 DATA "U6R4BR1BD1D1BD1D2BD1BL1L4U3R4"
6630 DATA "BU1U4BU1BR1R3F1BD4G1L3"
6640 DATA "U6R3F1F1D2G1G1L3"
6650 DATA "U6R5BD3BL4R3BD3BL4R5"
6660 DATA "U6R5BD3BL4R3"
6670 DATA "BU1U4BU1BR1R3F1BD2L2BR2D3L4"
6680 DATA "U6D3R5U3D6"
6690 DATA "BR2R2L1U6L1R2"
6700 DATA "BU1U1BF2L1R3BR1BU1U5BL1R2"
6710 DATA "U6D3R2E3G3F3"
6720 DATA "U6D6R5"
6730 DATA "U6F3D1U1E3D6"
6740 DATA "U6BD1R1BD1R1BD1R1BD1R1BD1R1BD1R1U6"
6750 DATA "BU1U4E1R3F1D4G1L3"
6760 DATA "U6R4F1D1G1L3"
6770 DATA "BU1U4BE1R3BF1D2G3L1BR2BU3F3"
6780 DATA "U6R4BF1D1BG1L3R1F3"
6790 DATA "BU1F1R3E1U1H1L3H1U1E1R3F1"
6800 DATA "BR2U6L3R6"
6810 DATA "BU1U5D5F1R3E1U5"
6820 DATA "BU3U3D3F3E3U3"
6830 DATA "U6D6E3U1D1F3U6"
6840 DATA "E6G3H3F6"
6850 DATA "BR3U3H3F3E3"
6860 DATA "R5L5E5U1L5"
6870 DATA "BR1R2L2U6R2"
6880 DATA "BU6F6"
6890 DATA "BR2R2U6L2"
6900 DATA "BR3U6F3H3G3"
6910 DATA "R6"
6920 DATA "BU6BR3D1G1"
6930 DATA "BR1H1E1R4D2L4R4U3H1L3"
6940 DATA "BR1R3E1U2H1L2G2D1U5"
6950 DATA "BR1R3E1G1L3H1U2E1R3"
6960 DATA "BR1R4U6D4H2L2G1D2"
6970 DATA "BR1R3L3H1U2E1R3F1D1L5"
6980 DATA "BR2U3L1R2L1U2E1R1F1"
6990 DATA "BR1R3E1U3L4G1F1R4"
7000 DATA "U6D4E2R2F1D3"
7010 DATA "BR1R2L1U4BU2L1BR1BD2L1"
7020 DATA "BU1F1R2E1U3BU2L1"
7030 DATA "BR1U6D4R2E2G2F2"
7040 DATA "BR2R2L1U6L1"
7050 DATA "U4R1R1F1D3U3E1R1F1D3"
7060 DATA "U4D2E2R2F1D3"
7070 DATA "BR1R3E1U2H1L3G1D2"
7080 DATA "U4R5F1G1L5"
7090 DATA "BR5U4L5G1F1R5"
7100 DATA "BR1U4D2E2R2"
7110 DATA "R4E1H1L3H1E1R4"
7120 DATA "BR2R1E1G1BL1BU1U4D1L1R2"
7130 DATA "BU1U3D3F1R2E3U1D4"
7140 DATA "BU3U1D1F3E3U1"
7150 DATA "BU1U3D3F1R1E1U1D1F1R1E1U3"
7160 DATA "E2R1E2G2L1H2F2R1F2"
7170 DATA "R2E3U1D1G1L2H2"
7180 DATA "R5L5BU1R1BR1BU1R1BU1BR1R1BU1L5"
7190 DATA "BR3R1L1H1U1H1L1R1E1U1E1R1"
7200 DATA "BR3U2BU2U2"
7210 DATA "BR2R1E1U1E1R1L1H1U1H1L1"
10000 FOR I = 1 TO N: TT(I) = TT(I) * 60: NEXT I: RETURN
***** End of Program Listing ******


Now there is also an Excel version of this program. Click here to download the Excel file or use the Visual Basic for Applications listing below.


***** Begin of Program Listing *****

'
'               NUMERICAL SMOOTHING AND DIFFERENTATION OF EMPIRICAL DATA
'
'
'     References:
'
'       - LE DUY, A. & ZAJIC, J.E.
'         A Geometrical Approach for Differentiation of Experimental Function at a Point
'         Biotechnology and Bioengineering, Vol. XV, 1973, 805-810
'
'       - AUBANEL, E.E. & OLDHAM, K.B.
'         Fourier Smoothing - Without the Fast Fourier Transform
'         Byte, February 1985, 207-218
'
'
'                 Antonio Augusto Gorni - www.gorni.eng.br
'
'                               Program Versions
'                - Applesoft Basic:               23.06.1988
'                - QuickBasic:                    14.10.1994
'                - Visual Basic for Applications: 14.01.2019
'

Option Explicit
Option Base 1

Const iLin As Integer = 6
Const Pi As Single = 3.141592654

Dim d1 As Integer
Dim d2 As Integer
Dim e As Integer                    ' Number of transform points to be kept
Dim i As Integer
Dim j As Integer
Dim j1 As Integer
Dim j2 As Integer
Dim k As Integer
Dim k1 As Integer
Dim l As Integer
Dim l1 As Integer
Dim n As Integer                    ' Number of data points.
Dim n2 As Integer
Dim p As Integer
Dim p1 As Integer

Dim b As Single
Dim c As Single
Dim EPct As Single                  ' Smoothing factor [%]
Dim f As Single
Dim f1 As Single
Dim f2 As Single
Dim g As Single
Dim m As Single
Dim q As Single
Dim s As Single
Dim s1 As Single
Dim s2 As Single
Dim u As Single
Dim t As Single
Dim t2 As Single
Dim x1 As Single
Dim x2 As Single

Dim Buf1 As String
Dim Buf2 As String
Dim DataCond As String

Dim Flag1 As Boolean

Dim Cr(5000) As Single              ' Internal array for derivative calculation.
Dim CX(5000) As Single              ' Dependent variable data for derivative calculation, original
Dim CXa(10000) As Single            ' Dependent variable data for derivative calculation, eventually smoothed
Dim DX(5000) As Single              ' Derivative of dependent variable.
Dim iSm(5000) As Single             ' i(k) transform for smoothing calculation.
Dim r(5000) As Single               ' r(k) transform for smoothing calculation.
Dim Tc(5000) As Single              ' Internal array for derivative calculation.
Dim Tt(5000) As Single              ' Independent variable data.
Dim uSm(5000) As Single             ' Internal array for smoothing calculation.
Dim v(5000) As Single               ' Internal array for smoothing calculation.
Dim Zb(5000) As Single              ' Internal array for derivative calculation.
Dim Zc(5000) As Single              ' Internal array for derivative calculation.
Dim Zi(5000) As Single              ' Internal array for derivative calculation.
Dim Zm(5000) As Single              ' Internal array for derivative calculation.
Dim Zn(5000) As Single              ' Internal array for derivative calculation.
Dim Zo(5000) As Single              ' Internal array for derivative calculation.

'
' Sub DataDiff: Calculation of the differentiation.
'

Sub DataDiff()

'
' Initialization and cleaning of previous results.
'

Application.ScreenUpdating = False

i = iLin
Do While Cells(i, "F") <> ""
   Cells(i, "F") = ""
   i = i + 1
Loop

'
' Input of data for differentiation.
'

Flag1 = True
Do While Flag1
   Buf1 = InputBox("Use smoothed data (Y/N)?", "Data")
   If Buf1 = "" Then End
   If Buf1 = "Y" Or Buf1 = "y" Then
      Buf2 = "D"
      Flag1 = False
   ElseIf Buf1 = "N" Or Buf1 = "n" Then
      Buf2 = "B"
      Flag1 = False
   End If
Loop

i = iLin
Do While Cells(i, "A") <> ""
   Tt(i - iLin + 1) = Cells(i, "A")
   CXa(i - iLin + 1) = Cells(i, Buf2)
   i = i + 1
Loop
n = i - iLin - 1

'
' Differentation of empirical data, eventually already smoothed.
'

For i = 1 To n - 2
   If CXa(i) = CXa(i + 1) Then CXa(i + 1) = CXa(i) * 0.9999
   If CXa(i) = CXa(i + 2) Then CXa(i + 2) = CXa(i + 1) * 0.9999
   If Tt(i) = Tt(i + 1) Then Tt(i + 1) = Tt(i) * 0.9999
   If Tt(i) = Tt(i + 2) Then Tt(i + 2) = Tt(i + 1) * 0.9999
Next i

If CXa(n) = CXa(n - 1) Then CXa(n) = CXa(n) * 0.9999
If Tt(n) = Tt(n - 1) Then Tt(n) = Tt(n) * 0.9999

i = 1
DX(i) = (CXa(i + 1) - CXa(i)) / (Tt(i + 1) - Tt(i))

Do While i - n <= -2
   Zb(i + 1) = (CXa(i + 1) - CXa(i)) / (Tt(i + 1) - Tt(i))
   Zc(i + 1) = (CXa(i + 2) - CXa(i + 1)) / (Tt(i + 2) - Tt(i + 1))
   If Abs(Zb(i + 1) - Zc(i + 1)) > 0.001 Then
      Zo(i + 1) = (Tt(i + 1) - Tt(i + 2)) / (CXa(i + 2) - CXa(i + 1))
      Zn(i + 1) = 0.5 * (CXa(i + 1) + CXa(i + 2)) - Zo(i + 1) * 0.5 * (Tt(i + 1) + Tt(i + 2))
      Zm(i + 1) = (Tt(i) - Tt(i + 1)) / (CXa(i + 1) - CXa(i))
      Zi(i + 1) = 0.5 * (CXa(i) + CXa(i + 1)) - Zm(i + 1) * 0.5 * (Tt(i) + Tt(i + 1))
      Tc(i + 1) = (Zn(i + 1) - Zi(i + 1)) / (Zm(i + 1) - Zo(i + 1))
      Cr(i + 1) = Zm(i + 1) * Tc(i + 1) + Zi(i + 1)
      DX(i + 1) = (Tc(i + 1) - Tt(i + 1)) / (CXa(i + 1) - Cr(i + 1))
                                      Else
      DX(i + 1) = (CXa(i + 1) - CXa(i)) / (Tt(i + 1) - Tt(i))
   End If
   If i < n Then DX(i + 1) = 0.5 * (Zb(i + 1) + Zc(i + 1))
   i = i + 1
Loop
   
'
' Output of differentation values.
'
   
For i = iLin To iLin + n - 1
   Cells(i, "F") = DX(i - iLin + 1)
Next

Application.ScreenUpdating = True

End Sub

'
' Subroutine SmoothPlain: Plain data smoothing.
'


Sub SmoothPlain()

DataCond = "Plain"
SmoothCentral (DataCond)

End Sub

'
' Subroutine SmoothRevMirror: Data smoothing but mirroring data starting at the final element of the set.
'


Sub SmoothRevMirror()

DataCond = "ReverseMirror"
SmoothCentral (DataCond)

End Sub

'
' Subroutine SmoothForwMirror: Data smoothing but mirroring data starting at the initial element of the set.
'


Sub SmoothForwMirror()

DataCond = "ForwardMirror"
SmoothCentral (DataCond)

End Sub

'
' Subroutine SmoothStrLine: Data smoothing but mirroring data starting at the initial element of the set.
'


Sub SmoothStrLine()

DataCond = "StraightLine"
SmoothCentral (DataCond)

End Sub



'
' Subroutine Smooth: Data smoothing, including options of data conditioning.
'

Sub SmoothCentral(DataCond As String)

'
' Initialization and cleaning of previous results.
'

Application.ScreenUpdating = False

i = iLin
Do While Cells(i, "F") <> ""
   Cells(i, "D") = ""
   Cells(i, "F") = ""
   i = i + 1
Loop

'
' Input of data for smoothing.
'

i = iLin
Do While Cells(i, "A") <> ""
   CX(i - iLin + 1) = Cells(i, "B")
   i = i + 1
Loop
n = i - iLin - 1

If DataCond = "Plain" Or DataCond = "ForwardMirror" Or DataCond = "StraightLine" Then
   For i = 1 To n
      CXa(i) = CX(i)
   Next i
End If

'
' Data processing to avoid smoothing distortions: reverse mirror.
'

If DataCond = "ReverseMirror" Then
   For i = n + 1 To 2 * n
      CXa(i) = CX(i - n)
   Next i
   For i = 1 To n
      CXa(i) = CX(n - i + 1)
   Next i
   n = 2 * n
End If

'
' Data processing to avoid smoothing distortions: forward mirror.
'

If DataCond = "ForwardMirror" Then
   For i = n + 1 To 2 * n
      CXa(i) = CX(2 * n - i + 1)
   Next i
   n = 2 * n
End If

'
' Data processing to avoid smoothing distortions: approximation using straight line.
'

If DataCond = "StraightLine" Then
   d1 = Int(n / 10)
   d2 = d1
   d1 = InputBox("No. of Reference Points, Lowest Extremity?" & Chr(13) & "(Default =" & d1 & ")", "Data")
   d2 = InputBox("No. of Reference Points, Highest Extremity?" & Chr(13) & "(Default =" & d2 & ")", "Data")
End If

'
' User must define smoothing factor.
'

q = 0
n2 = Int((n + 1) / 2 + 1)
EPct = InputBox("Smoothing Factor [%]?", "Data")
e = Int((100 - EPct) * (n + 1) / 200)
If e > Int((n + 1) / 2) Then
   MsgBox "Insuficient Smoothing Factor! Increase it!", 0, "Warning!"
   Exit Sub
End If
If e <= 1 Then
   MsgBox "Excessive Smoothing Factor! Decrease it!", 0, "Warning!"
   Exit Sub
End If

'
' Smoothing calculation start.
'

If e > q Then
   If q = 0 Then
      
'
' Optional straight-line procedure for eliminating end effect was chosen?
'
      
      If DataCond = "StraightLine" Then
         
'
' Yes, perform straight-line procedure.
'
         
         s1 = 0
         s2 = 0
         For j = 0 To d1 - 1
            s1 = s1 + CXa(j + 1)
         Next j
         For j = 0 To d2 - j
            s2 = s2 + CXa(n - j)
         Next j
         x1 = s1 / d1
         x2 = s2 / d2
         m = (x2 - x1) / (n - (d1 + d2) / 2)
         b = (x1 + x2) / 2 - m * n / 2
                                    
'
' No, straight-line procedure was not chosen.
'
                                    
                                    Else
         m = 0
         b = 0
      End If

'
' Calculate 'r(0)' transform.
'

      g = 0
      For j = 0 To n - 1
         CXa(j + 1) = CXa(j + 1) - m * j - b
         g = g + CXa(j + 1)
      Next j
      r(1) = g / n
      q = 1
   End If
   
'
' Calculate 'r(k)' transforms.
'
   
   j2 = Int((n - 1) / 2)
   p1 = Int(Log(2 * j2 - 1) / Log(2))
   For k = q To e - 1
      j1 = j2
      s = Pi * k * 2 / n
      c = Cos(s)
      s = Sin(s)
      For j = 1 To j1
         l = 2 * j - 1
         uSm(j + 1) = CXa(l + 1) * c + CXa(l + 2)
         v(j + 1) = CXa(l + 1) * s
      Next j
      s = 2 * s * c
      c = 2 * c * c - 1
      For p = 1 To p1
         uSm(j1 + 2) = 0
         v(j1 + 2) = 0
         j1 = Int((j1 + 1) / 2)
         For j = 1 To j1
            l = 2 * j - 1
            u = uSm(l + 1) * c - v(l + 1) * s + uSm(l + 2)
            v(j + 1) = uSm(l + 1) * s + v(l + 1) * c + v(l + 2)
            uSm(j + 1) = u
         Next j
         s = 2 * s * c
         c = 2 * c * c - 1
      Next p
      r(k + 1) = (CXa(1) + (uSm(2) * c + v(2) * s)) / n
   Next k
   
'
' Calculate 'iSm(k)' Transforms
'
   
   For k = q To e - 1
      j1 = j2
      s = 2 * Pi * k / n
      c = Cos(s)
      s = Sin(s)
      For j = 1 To j1
         l = 2 * j - 1
         uSm(j + 1) = -(CXa(l + 1) * s)
         v(j + 1) = CXa(l + 1) * c + CXa(l + 2)
      Next j
      s = 2 * s * c
      c = 2 * c * c - 1
      For p = 1 To p1
         uSm(j1 + 2) = 0
         v(j1 + 2) = 0
         j1 = Int((j1 + 1) / 2)
         For j = 1 To j1
            l = 2 * j - 1
            u = uSm(l + 1) * c - v(l + 1) * s + uSm(l + 2)
            v(j + 1) = uSm(l + 1) * s + v(l + 1) * c + v(l + 2)
            uSm(j + 1) = u
         Next j
         s = 2 * s * c
         c = 2 * c * c - 1
      Next p
      iSm(k + 1) = -((uSm(2) * c + v(2) * s) / n)
   Next k
   If e > q Then q = e
End If

'
' Working on inverse transform.
'
' Calculate 'x1(0)'.
'

f1 = 0
f2 = 0
For k = 1 To e - 1
   t = r(k + 1)
   f1 = f1 + t
   f2 = f2 + CSng(k) * CSng(k) * t
Next k
CXa(1) = r(1) + 2 * (f1 - f2 * (1 / e / e))

'
' Next instruction is to be performed if straight-line procedure was chosen.
'

If DataCond = "StraightLine" Then CXa(1) = CXa(1) + b

p1 = Int(Log(2 * e - 3) / Log(2))
For j = 1 To n - 1
   t2 = CSng(e) * CSng(e)
   For k = 1 To e - 1
      f = 1 - CSng(k) * CSng(k) / t2
      uSm(k + 1) = r(k + 1) * f
      v(k + 1) = -(iSm(k + 1) * f)
   Next k
   k1 = e - 1
   s = 2 * Pi * j / n
   c = Cos(s)
   s = Sin(s)
   For p = 1 To p1
      uSm(k1 + 2) = 0
      v(k1 + 2) = 0
      k1 = Int((k1 + 1) / 2)
      For k = 1 To k1
         l = 2 * k - 1
         u = uSm(l + 1) * c - v(l + 1) * s + uSm(l + 2)
         v(k + 1) = uSm(l + 1) * s + v(l + 1) * c + v(l + 2)
         uSm(k + 1) = u
      Next k
      s = 2 * s * c
      c = 2 * c * c - 1
   Next p
   CXa(j + 1) = r(1) + 2 * (uSm(2) * c + v(2) * s)
   
'
' Next instruction is to be performed if straight-line procedure was chosen.
'
   
   If DataCond = "StraightLine" Then CXa(j + 1) = CXa(j + 1) + m * j + b
Next j
 
'
' Smoothing calculating end.
'
' Is data "demirroring" needed?
'
 
If DataCond = "ReverseMirror" Or DataCond = "ForwardMirror" Then
   
'
' Yes, proceed.
'
   
   n = n / 2
   If DataCond = "ReverseMirror" Then
      For i = 1 To n
         CXa(i) = CXa(i + n)
      Next i
   End If
End If

'
' Output of smoothed data.
'

For i = 1 To n
   Cells(i + iLin - 1, "D") = CXa(i)
Next i

Application.ScreenUpdating = True

End Sub

'
' Sub Clean: Clears all data and results.
'

Sub Clean()

Application.ScreenUpdating = False

i = iLin
Do While Cells(i, "A") <> ""
   Rows(i).Select
   Selection.ClearContents
   i = i + 1
Loop

Application.ScreenUpdating = True

End Sub

***** End of Program Listing ******


Return to the Software Menu.

Originally Posted: 25 June 1996 - Last Update: 19 January 2019
© Antonio Augusto Gorni