summaryrefslogtreecommitdiff
path: root/romunit.pas
diff options
context:
space:
mode:
Diffstat (limited to 'romunit.pas')
-rw-r--r--romunit.pas202
1 files changed, 1 insertions, 201 deletions
diff --git a/romunit.pas b/romunit.pas
index c1c0a7c..761bd68 100644
--- a/romunit.pas
+++ b/romunit.pas
@@ -10,8 +10,7 @@ uses
procedure Fehler(s: string);
procedure readRawInputs(nam: string; out d1,d2: tExtPointArray; var absShift: extended);
procedure readTextInput(nam: string; out dat: tExtPointArray);
-procedure writeOutput(nam: string; const dat: tExtPointArray); overload;
-procedure writeOutput(nam: string; const dat,datArg: tExtPointArray); overload;
+procedure writeOutput(nam: string; const dat: tExtPointArray);
procedure berechneTrajektorie(const inPuls,outPuls: tExtPointArray; out surTraj: tExtPointArray; absShift: extended);
procedure smoothen(var dat: tExtPointArray; width: longint);
procedure sort(var dat: tExtPointArray); overload;
@@ -28,11 +27,9 @@ procedure removeLinearOffset(var dat: tExtPointArray);
procedure flip(var dat: tExtPointArray);
procedure vereineExtrema(const dat1,dat2: tExtPointArray; var el1,el2: tLongintArray; toleranz: extended);
procedure uniq(var dat: tExtPointArray; streng: boolean);
-procedure fft(var dat: tExtPointArray; out datArg: tExtPointArray);
procedure interpoliere(var dat: tExtPointArray);
procedure normiere(var dat: tExtPointArray);
procedure berechneRefPuls(inPuls,surTraj: tExtPointArray; betaGlaette: longint; betaBound: extended; out cRefPuls: tExtPointArray);
-procedure fenstern(var dat: tExtPointArray; breite: extended; alteMitte: boolean = false);
type
tSortThread = class(tThread)
@@ -306,34 +303,6 @@ begin
writeln(stderr,' fertig');
end;
-procedure writeOutput(nam: string; const dat,datArg: tExtPointArray); overload;
-var
- f: textfile;
- i: longint;
-begin
- if length(dat)=0 then
- exit;
- if length(dat)<>length(datArg) then
- Fehler('dat und datArg haben unterschiedlich viele Werte!');
- write(stderr,'Datei '''+nam+''' schreiben ');
- assignfile(f,nam);
- rewrite(f);
- for i:=0 to length(dat)-1 do begin
- if dat[i]['x']<>datArg[i]['x'] then
- Fehler('dat und datArg haben unterschiedliche x-Werte an Position '+inttostr(i)+' ('+floattostr(dat[i]['x'])+' vs. '+floattostr(datArg[i]['x'])+')!');
- writeln(
- f,
- floattostr(dat[i]['x'])+#9+
- floattostr(dat[i]['y'])+#9+
- floattostr(datArg[i]['y'])
- );
- if i and 65535 = 0 then
- write(stderr,'.');
- end;
- closefile(f);
- writeln(stderr,' fertig');
-end;
-
procedure berechneTrajektorie(const inPuls,outPuls: tExtPointArray; out surTraj: tExtPointArray; absShift: extended);
var
i,j: longint;
@@ -811,119 +780,6 @@ begin
setlength(dat,j);
end;
-procedure fft(var dat: tExtPointArray; out datArg: tExtPointArray);
-var
- i,j,k,n,dist,absch,wnum,wstep,haL: longint;
- in0,out0: boolean;
- wRe,wIm: tExtendedArray;
- t1,t2,vorher,nachher,fstep,pvFehler: extended;
- umsortierung: tLongintArray;
-begin
- write(stderr,'FFT ('+inttostr(length(dat))+') ... ');
- fstep:=dat[1]['x']-dat[0]['x'];
- j:=length(dat);
- setlength(dat,2*round(power(2,ceil(ln(length(dat))/ln(2)))));
- k:=(length(dat)-j) div 2;
- for i:=j-1 downto 0 do
- dat[i+k]:=dat[i];
- for i:=0 to k-1 do // weich auf 0 abklingen lassen
- dat[i]['y']:=dat[2*k-1-i]['y'] * sqr(sin(pi/2*i/k));
- for i:=j+k to length(dat)-1 do // dito
- dat[i]['y']:=dat[2*j+2*k-1-i]['y'] * sqr(sin(pi/2*(length(dat)-1-i)/(length(dat)-j-k)));
- for i:=0 to length(dat)-1 do begin
- dat[i]['x']:=dat[i]['y'];
- dat[i]['y']:=0;
- end;
-
- haL:=length(dat) div 2;
- n:=round(ln(2*hal)/ln(2));
-
- fstep:=1/(2*haL)/fstep;
-
- vorher:=0;
- for i:=0 to 2*haL-1 do
- vorher:=vorher + sqr(dat[i]['x']);
-
- setlength(umsortierung,2*haL);
- for j:=0 to 2*haL-1 do begin
- k:=0;
- for i:=0 to n-1 do // Bits umkehren
- k:=2*k + byte(odd(j shr i));
- umsortierung[j]:=k;
- end;
- setlength(wRe,haL);
- setlength(wIm,haL);
- for i:=0 to haL-1 do begin
- wRe[i]:=cos(pi*i/haL);
- wIm[i]:=sin(pi*i/haL);
- end;
-
- in0:=true;
- out0:=true;
- for i:=0 to 2*haL-1 do begin
- if umsortierung[i]>i then begin
- t1:=dat[i]['x'];
- dat[i]['x']:=dat[umsortierung[i]]['x'];
- dat[umsortierung[i]]['x']:=t1;
- end;
- in0:=in0 and (dat[i]['x']=0);
- end;
-
- dist:=1;
- wstep:=haL;
- while dist<2*haL do begin
- absch:=0;
- while absch<2*haL do begin
- wnum:=0;
- for i:=0 to dist-1 do begin
- // x_links: [absch+j]
- // x_rechts: [absch+j+dist]
- t1:=wRe[wnum]*dat[absch+i+dist]['x'] - wIm[wnum]*dat[absch+i+dist]['y'];
- t2:=wRe[wnum]*dat[absch+i+dist]['y'] + wIm[wnum]*dat[absch+i+dist]['x'];
-
- dat[absch+i+dist]['x']:=dat[absch+i]['x']-t1;
- dat[absch+i+dist]['y']:=dat[absch+i]['y']-t2;
- dat[absch+i]['x']:=dat[absch+i]['x']+t1;
- dat[absch+i]['y']:=dat[absch+i]['y']+t2;
- wnum:=wnum+wstep;
- end;
- absch:=absch+2*dist;
- end;
- wstep:=wstep div 2;
- dist:=dist*2;
- end;
-
- setlength(datArg,length(dat));
-
- t1:=0;
- for i:=0 to 2*haL-1 do begin
- t2:=arctan2(dat[i]['y'],dat[i]['x']);
- t2:=t2-2*pi*round((t2-t1)/2/pi);
- datArg[i]['y']:=t2;
- t1:=t2;
- dat[i]['y']:=(sqr(dat[i]['x'])+sqr(dat[i]['y']))/(2*haL);
- dat[i]['x']:=fStep*i;
- datArg[i]['x']:=dat[i]['x'];
- end;
- for i:=0 to 2*haL-1 do
- out0:=out0 and (dat[i]['y']=0);
-
- nachher:=0;
- for i:=0 to 2*haL-1 do
- nachher:=nachher + dat[i]['y'];
-
- if (nachher=0) and (vorher=0) then
- pvFehler:=0
- else
- pvFehler:=abs(nachher-vorher)/(nachher+vorher);
-
- writeln(stderr,' fertig (Parseval-Fehler = '+floattostr(pvFehler)+' ... nämlich von '+floattostr(vorher)+' zu '+floattostr(nachher)+')');
- if in0 then
- writeln(stderr,'Nur Nullen im Input der FFT!');
- if out0 then
- writeln(stderr,'Nur Nullen im Output der FFT!');
-end;
-
procedure interpoliere(var dat: tExtPointArray);
var
i,j: longint;
@@ -1121,61 +977,5 @@ begin
fertig:=true;
end;
-var
- fensterMitte: longint;
-
-procedure fenstern(var dat: tExtPointArray; breite: extended; alteMitte: boolean = false);
-var
- i,j,len: longint;
- extrema: tLongintArray;
- pps: tExtendedArray;
- mpp: extended;
-begin
- if not alteMitte then begin
- setlength(extrema,0);
- len:=0;
- for i:=1 to length(dat)-2 do
- if (dat[i]['y']>dat[i-1]['y']) xor
- (dat[i+1]['y']>dat[i]['y']) then begin
- if length(extrema)<=len then
- setlength(extrema,len+65536);
- extrema[len]:=i;
- inc(len);
- end;
- setlength(extrema,len);
- setlength(pps,len);
- for i:=0 to length(extrema)-1 do begin
- pps[i]:=-1;
- for j:=0 to length(extrema)-1 do
- if (abs(dat[extrema[j]]['x']-dat[extrema[i]]['x'])<2) and
- (abs(dat[extrema[j]]['y']-dat[extrema[i]]['y'])>pps[i]) then
- pps[i]:=abs(dat[extrema[j]]['y']-dat[extrema[i]]['y']);
- end;
- mpp:=pps[0];
- for i:=0 to length(pps)-1 do
- if pps[i]>mpp then
- mpp:=pps[i];
- i:=0;
- while (i<length(pps)) and
- (pps[i]<mpp/2) do
- inc(i);
- j:=length(pps)-1;
- while (j>=0) and
- (pps[j]<mpp/2) do
- dec(j);
- mpp:=(dat[extrema[i]]['x']+dat[extrema[j]]['x'])/2;
- fensterMitte:=0;
- while (fensterMitte<length(dat)) and
- (dat[fensterMitte]['x']<mpp) do
- inc(fensterMitte);
- end;
-
- for i:=0 to length(dat)-1 do
- if abs(dat[i]['x']-dat[fensterMitte]['x']) >= breite then
- dat[i]['y']:=0
- else
- dat[i]['y']:=dat[i]['y'] * (1+cos((dat[i]['x']-dat[fensterMitte]['x'])/breite*pi))/2;
-end;
-
end.