From 8b3a90a0c6cc926428844f4474db4fa674d1933b Mon Sep 17 00:00:00 2001 From: Erich Eckner Date: Tue, 19 Jun 2018 11:06:46 +0200 Subject: tExtPoint is now an array instead of a record --- ROM.lpi | 5 +- ROM.lpr | 16 ++-- ROM.lps | 112 +++++++++++++-------------- romunit.pas | 245 ++++++++++++++++++++++++++++++------------------------------ 4 files changed, 188 insertions(+), 190 deletions(-) diff --git a/ROM.lpi b/ROM.lpi index 07e671c..0e77c85 100644 --- a/ROM.lpi +++ b/ROM.lpi @@ -1,7 +1,7 @@ - + @@ -16,9 +16,6 @@ - - - diff --git a/ROM.lpr b/ROM.lpr index 898a35b..941a700 100644 --- a/ROM.lpr +++ b/ROM.lpr @@ -10,7 +10,7 @@ uses {$ENDIF}{$ENDIF} Classes { you can add units after this }, - SysUtils,ROMunit, matheunit, Math, systemunit; + SysUtils,ROMunit, matheunit, Math, systemunit, lowlevelunit; var inPulsO,inPuls,refPulsO,refPuls, @@ -277,14 +277,14 @@ begin fft(refPuls); fft(surTraj); fft(surVel); - inPuls[0].y:=0; - refPuls[0].y:=0; - surTraj[0].y:=0; + inPuls[0]['y']:=0; + refPuls[0]['y']:=0; + surTraj[0]['y']:=0; if wmax<0 then begin - cut(surTraj,min(refPuls[length(refPuls)-1].x,inPuls[length(inPuls)-1].x)/2); - cut(inPuls,surTraj[length(surTraj)-1].x); - cut(refPuls,surTraj[length(surTraj)-1].x); - cut(surVel,surTraj[length(surTraj)-1].x); + cut(surTraj,min(refPuls[length(refPuls)-1]['x'],inPuls[length(inPuls)-1]['x'])/2); + cut(inPuls,surTraj[length(surTraj)-1]['x']); + cut(refPuls,surTraj[length(surTraj)-1]['x']); + cut(surVel,surTraj[length(surTraj)-1]['x']); end else begin cut(surTraj,wmax); diff --git a/ROM.lps b/ROM.lps index 77633ee..8e5b84a 100644 --- a/ROM.lps +++ b/ROM.lps @@ -1,24 +1,26 @@ - + - - - - + + + + - - + + + + @@ -26,14 +28,12 @@ - + - - - + @@ -45,128 +45,130 @@ - - + - - + + - - + + - + - + - + - - + + - - + + - + + - + + - + - - + + - + - + + - - + + - - + + - - + + - - + + - + - + - + - - + + - + - + - + - - + + - + - - + + - - + + - + diff --git a/romunit.pas b/romunit.pas index 9cf905e..6bce2e0 100644 --- a/romunit.pas +++ b/romunit.pas @@ -5,7 +5,7 @@ unit romunit; interface uses - Classes, SysUtils, Math, matheunit; + Classes, SysUtils, Math, matheunit, lowlevelunit; procedure Fehler(s: string); procedure readRawInputs(nam: string; out d1,d2: tExtPointArray; var absShift: extended); @@ -193,13 +193,13 @@ begin if i=start then cells_left:=cells_left-fl; blockread(f,buff[0],sizeof(single)*steps); // fp for k:=0 to steps-1 do begin - d1[k+length(d1)-steps].x:=(i-1+k/steps)*factor; - d1[k+length(d1)-steps].y:=buff[k]/factor; + d1[k+length(d1)-steps]['x']:=(i-1+k/steps)*factor; + d1[k+length(d1)-steps]['y']:=buff[k]/factor; end; blockread(f,buff[0],sizeof(single)*steps); // fm for k:=0 to steps-1 do begin - d2[k+length(d2)-steps].x:=(i-1+k/steps)*factor; - d2[k+length(d2)-steps].y:=buff[k]/factor; + d2[k+length(d2)-steps]['x']:=(i-1+k/steps)*factor; + d2[k+length(d2)-steps]['y']:=buff[k]/factor; end; end; for k:=0 to 7+2*byte(j<>0) do // (fp,fm,)gp,gm,ex,de,di,jx,jy,jz @@ -244,10 +244,10 @@ begin if pos('#',s)>0 then delete(s,pos('#',s),length(s)); while (length(s)>0) and (s[1]=' ') do delete(s,1,1); if s='' then continue; - dat[i].x:=strtofloat(copy(s,1,pos(' ',s)-1)); + dat[i]['x']:=strtofloat(copy(s,1,pos(' ',s)-1)); delete(s,1,pos(' ',s)); while (length(s)>0) and (s[1]=' ') do delete(s,1,1); - dat[i].y:=strtofloat(copy(s,1,pos(' ',s+' ')-1)); + dat[i]['y']:=strtofloat(copy(s,1,pos(' ',s+' ')-1)); if i and 1023 = 0 then write(stderr,'.'); end; closefile(f); @@ -268,8 +268,8 @@ begin ml:=0; mr:=0; for i:=0 to length(dat)-1 do begin - strings[i,false]:=floattostr(dat[i].x); - strings[i,true]:=floattostr(dat[i].y); + strings[i,false]:=floattostr(dat[i]['x']); + strings[i,true]:=floattostr(dat[i]['y']); ml:=max(ml,length(strings[i,false])); mr:=max(mr,pos('.',strings[i,true]+'.')); if i and 65535 = 0 then write(stderr,'.'); @@ -293,10 +293,10 @@ var i,j: longint; b1,b2: boolean; begin for i:=0 to length(inPuls)-2 do - if inPuls[i+1].x<=inPuls[i].x then + if inPuls[i+1]['x']<=inPuls[i]['x'] then writeln('Der Inpuls ist nicht sortiert!'); for i:=0 to length(outPuls)-2 do - if outPuls[i+1].x<=outPuls[i].x then + if outPuls[i+1]['x']<=outPuls[i]['x'] then writeln('Der Outpuls ist nicht sortiert!'); // Extrema identifizieren @@ -438,7 +438,7 @@ begin for i:=0 to length(apos)-1 do begin if apos[i]>epos[i] then continue; fertig:=false; - if (j=-1) or (dat[apos[i]].xdat[i-1].y) and (dat[i].y>=dat[i+1].y) and (dat[i].y>0)) or - ((not maxima) and (dat[i].ydat[i-1]['y']) and (dat[i]['y']>=dat[i+1]['y']) and (dat[i]['y']>0)) or + ((not maxima) and (dat[i]['y']dat[i-1].y) and (dat[i].y>=dat[i+1].y) and (dat[i].y>0)) or - ((not maxima) and (dat[i].ydat[i-1]['y']) and (dat[i]['y']>=dat[i+1]['y']) and (dat[i]['y']>0)) or + ((not maxima) and (dat[i]['y']mx then begin - mx:=abs(dat[extrema[i]].y);; + if abs(dat[extrema[i]]['y'])>mx then begin + mx:=abs(dat[extrema[i]]['y']);; mitte:=i; end; - mx:=dat[extrema[mitte]].x; + mx:=dat[extrema[mitte]]['x']; - for i:=round((dat[extrema[0]].x - mx)/sollAbst) to round((dat[extrema[length(extrema)-1]].x - mx)/sollAbst) do begin - dist:=abs(dat[extrema[0]].x - i*sollAbst - mx); + for i:=round((dat[extrema[0]]['x'] - mx)/sollAbst) to round((dat[extrema[length(extrema)-1]]['x'] - mx)/sollAbst) do begin + dist:=abs(dat[extrema[0]]['x'] - i*sollAbst - mx); k:=0; for j:=1 to length(behalten)-1 do - if abs(dat[extrema[j]].x - i*sollAbst - mx) < dist then begin - dist:=abs(dat[extrema[j]].x - i*sollAbst - mx); + if abs(dat[extrema[j]]['x'] - i*sollAbst - mx) < dist then begin + dist:=abs(dat[extrema[j]]['x'] - i*sollAbst - mx); k:=j; end; if 2*dist=length(minima)) or ((idat[maxima[k]].y) then + if (idat[maxima[k]]['y']) then k:=i; inc(i); until (i>=length(maxima)) or ((j=minima[j])); @@ -584,22 +584,22 @@ var iMax,oMax,i: longint; begin iMax:=0; for i:=1 to length(inPuls)-1 do - if inPuls[i].y>inPuls[iMax].y then iMax:=i; + if inPuls[i]['y']>inPuls[iMax]['y'] then iMax:=i; oMax:=0; for i:=1 to length(outPuls)-1 do - if outPuls[i].y>outPuls[oMax].y then oMax:=i; - if absShift<-0.9e9 then absShift:=outPuls[oMax].x-inPuls[iMax].x - else absShift:=(outPuls[oMax].x-inPuls[iMax].x)-round((outPuls[oMax].x-inPuls[iMax].x)-absShift); + if outPuls[i]['y']>outPuls[oMax]['y'] then oMax:=i; + if absShift<-0.9e9 then absShift:=outPuls[oMax]['x']-inPuls[iMax]['x'] + else absShift:=(outPuls[oMax]['x']-inPuls[iMax]['x'])-round((outPuls[oMax]['x']-inPuls[iMax]['x'])-absShift); for i:=0 to length(outPuls)-1 do - outPuls[i].x:=outPuls[i].x-absShift; + outPuls[i]['x']:=outPuls[i]['x']-absShift; oMax:=0; - while (oMaxbound then - outdat[i].y:=nan; + outdat[i]['x']:=(indat[i+dist]['x']+indat[i]['x'])/2; +// outdat[i]['y']:=min(max((indat[i+dist]['y']-indat[i]['y'])/(indat[i+dist]['x']-indat[i]['x']),-bound),bound); + outdat[i]['y']:=(indat[i+dist]['y']-indat[i]['y'])/(indat[i+dist]['x']-indat[i]['x']); + if outdat[i]['y']>bound then + outdat[i]['y']:=bound; + if outdat[i]['y']<-bound then + outdat[i]['y']:=-bound; end; end; procedure removeLinearOffset(var dat: tExtPointArray); -var i: longint; - dx,dy: extended; +var i: longint; + d: tExtPoint; begin - dx:=0; - dy:=0; - for i:=0 to 99 do begin - dy:=dy + dat[length(dat)-1-i].y - dat[i].y; - dx:=dx + dat[length(dat)-1-i].x - dat[i].x; - end; + d:=extPoint(0,0); + for i:=0 to 99 do + d:=d + dat[length(dat)-1-i] - dat[i]; for i:=0 to length(dat)-1 do - dat[i].y:=dat[i].y - dy/dx*(dat[i].x-dat[0].x); - dy:=0; + dat[i]['y']:=dat[i]['y'] - d['y']/d['x']*(dat[i]['x']-dat[0]['x']); + d['y']:=0; for i:=0 to 99 do - dy:=dy + dat[length(dat)-1-i].y + dat[i].y; - dy:=dy/200; + d['y']:=d['y'] + dat[length(dat)-1-i]['y'] + dat[i]['y']; + d['y']:=d['y']/200; for i:=0 to length(dat)-1 do - dat[i].y:=dat[i].y - dy; + dat[i]['y']:=dat[i]['y'] - d['y']; end; procedure flip(var dat: tExtPointArray); var i: longint; begin for i:=0 to length(dat)-1 do - dat[i].y:=-dat[i].y; + dat[i]['y']:=-dat[i]['y']; end; procedure vereineExtrema(const dat1,dat2: tExtPointArray; var el1,el2: tLongintArray; toleranz: extended); @@ -680,8 +679,8 @@ begin for i:=0 to length(el1)-1 do begin k:=-1; for j:=0 to length(el2)-1 do - if (abs(dat1[el1[i]].x-dat2[el2[j]].x)<=toleranz) and - ((k=-1) or (abs(dat1[el1[i]].x-dat2[el2[j]].x)-1) and not behalten[true,k] then begin behalten[false,i]:=true; @@ -715,14 +714,14 @@ begin k:=0; tmp:=0; for i:=0 to length(dat)-1 do begin - if (i0) then begin - dat[j].y:=(tmp+dat[i].y)/(k+1); - dat[j].x:=dat[i].x; + dat[j]['y']:=(tmp+dat[i]['y'])/(k+1); + dat[j]['x']:=dat[i]['x']; k:=0; tmp:=0; end @@ -743,19 +742,19 @@ var i,j,k,n,dist,absch,wnum,wstep,haL: longint; umsortierung: tLongintArray; begin write(stderr,'FFT ... '); - fstep:=dat[1].x-dat[0].x; + 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[k].y * sqr(sin(pi/2*i/k)); + dat[i]['y']:=dat[k]['y'] * sqr(sin(pi/2*i/k)); for i:=j+k to length(dat)-1 do // dito - dat[i].y:=dat[j+k-1].y * sqr(sin(pi/2*(length(dat)-1-i)/(length(dat)-j-k))); + dat[i]['y']:=dat[j+k-1]['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; + dat[i]['x']:=dat[i]['y']; + dat[i]['y']:=0; end; haL:=length(dat) div 2; @@ -765,7 +764,7 @@ begin vorher:=0; for i:=0 to 2*haL-1 do - vorher:=vorher + sqr(dat[i].x); + vorher:=vorher + sqr(dat[i]['x']); setlength(umsortierung,2*haL); for j:=0 to 2*haL-1 do begin @@ -785,11 +784,11 @@ begin 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; + t1:=dat[i]['x']; + dat[i]['x']:=dat[umsortierung[i]]['x']; + dat[umsortierung[i]]['x']:=t1; end; - in0:=in0 and (dat[i].x=0); + in0:=in0 and (dat[i]['x']=0); end; dist:=1; @@ -801,13 +800,13 @@ begin 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; + 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; + 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; @@ -817,15 +816,15 @@ begin end; for i:=0 to 2*haL-1 do begin - dat[i].y:=(sqr(dat[i].x)+sqr(dat[i].y))/(2*haL); - dat[i].x:=fStep*i; + dat[i]['y']:=(sqr(dat[i]['x'])+sqr(dat[i]['y']))/(2*haL); + dat[i]['x']:=fStep*i; end; - for j:=0 to 2*haL-1 do - out0:=out0 and (dat[i].y=0); + 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; + nachher:=nachher + dat[i]['y']; if (nachher=0) and (vorher=0) then pvFehler:=0 else pvFehler:=abs(nachher-vorher)/(nachher+vorher); @@ -840,25 +839,25 @@ var i,j: longint; tmp: extended; tdat: tExtPointArray; begin - tmp:=dat[1].x-dat[0].x; + tmp:=dat[1]['x']-dat[0]['x']; for i:=2 to length(dat)-1 do - tmp:=min(tmp,dat[i].x-dat[i-1].x); + tmp:=min(tmp,dat[i]['x']-dat[i-1]['x']); if tmp<=0 then Fehler('Die Daten müssen sortiert sein und dürfen keine doppelten x-Werte enthalten! ('+floattostr(tmp)+')'); setlength(tdat,length(dat)); for i:=0 to length(dat)-1 do tdat[i]:=dat[i]; - setlength(dat,max(2*length(tdat),round(min(power(length(tdat),1.3),(tdat[length(tdat)-1].x-tdat[0].x)/tmp+1)))); + setlength(dat,max(2*length(tdat),round(min(power(length(tdat),1.3),(tdat[length(tdat)-1]['x']-tdat[0]['x'])/tmp+1)))); j:=0; for i:=0 to length(dat)-1 do begin - dat[i].x:=tdat[0].x+(tdat[length(tdat)-1].x-tdat[0].x)*i/(length(dat)-1); - while (jsurTraj[iSur].x-surTraj[iSur].y then + if inPuls[0]['x']>surTraj[iSur]['x']-surTraj[iSur]['y'] then continue; - while (iIn=length(inPuls)-1 then break; if anz>=length(cRefPuls) then setlength(cRefPuls,anz+32768); - beta:=min(max((surTraj[iSur+betaGlaette].y-surTraj[iSur].y)/(surTraj[iSur+betaGlaette].x-surTraj[iSur].x),-betaBound),betaBound); + beta:=min(max((surTraj[iSur+betaGlaette]['y']-surTraj[iSur]['y'])/(surTraj[iSur+betaGlaette]['x']-surTraj[iSur]['x']),-betaBound),betaBound); // beta:=0; - frac:=(surTraj[iSur].x-surTraj[iSur].y - inPuls[iIn].x)/(inPuls[iIn+1].x - inPuls[iIn].x); - cRefPuls[anz].x:=surTraj[iSur].x+surTraj[iSur].y; - cRefPuls[anz].y:= - (inPuls[iIn].y * (1-frac) + - inPuls[iIn+1].y * frac) * + frac:=(surTraj[iSur]['x']-surTraj[iSur]['y'] - inPuls[iIn]['x'])/(inPuls[iIn+1]['x'] - inPuls[iIn]['x']); + cRefPuls[anz]['x']:=surTraj[iSur]['x']+surTraj[iSur]['y']; + cRefPuls[anz]['y']:= + (inPuls[iIn]['y'] * (1-frac) + + inPuls[iIn+1]['y'] * frac) * (1-beta)/(1+beta); inc(anz); end; @@ -965,34 +964,34 @@ begin t1:=0; t2:=0; for i:=0 to length(_dats[b])-1 do begin - t1:=max(t1,_dats[b,i].y); - t2:=min(t2,_dats[b,i].y); + t1:=max(t1,_dats[b,i]['y']); + t2:=min(t2,_dats[b,i]['y']); end; for i:=0 to length(_dats[b])-1 do - _dats[b,i].y:=(_dats[b,i].y-t2)/(t1-t2); + _dats[b,i]['y']:=(_dats[b,i]['y']-t2)/(t1-t2); end; for b:=false to true do begin offset:=byte(b)*(_stops[false]-_starts[false]+1); for i:=0 to _stops[b]-_starts[b] do begin - _ergebnis[offset + i].x:=_dats[b,i].x; - if (i=0) then bestPos:=_dats[not b,0].x + _ergebnis[offset + i]['x']:=_dats[b,i]['x']; + if (i=0) then bestPos:=_dats[not b,0]['x'] else begin - t1:=_dats[b,i].y; - j:=0; // floor(_ergebnis[offset + i - 1].y); <- das war Murks + t1:=_dats[b,i]['y']; + j:=0; // floor(_ergebnis[offset + i - 1]['y']); <- das war Murks (* while (j<_stops[not b]-_starts[not b]) and - (_dats[not b,j].x<_ergebnis[offset + i - 1].y) do + (_dats[not b,j]['x']<_ergebnis[offset + i - 1]['y']) do inc(j);*) bestPos:=-1; bestErr:=-1; repeat while (j<_stops[not b]-_starts[not b]) and - ((_dats[not b,j].y - t1) * (_dats[not b,j + 1].y - t1) > 0) do + ((_dats[not b,j]['y'] - t1) * (_dats[not b,j + 1]['y'] - t1) > 0) do inc(j); - t2:=_dats[not b,j].x; + t2:=_dats[not b,j]['x']; if j<_stops[not b]-_starts[not b] then // x1 = dx * y1 / (y1 - y2) ... Berechnung des gebrochenen Anteils - t2:=t2 + (_dats[not b,j+1].x-_dats[not b,j].x) * (_dats[not b,j].y - t1) / (_dats[not b,j].y - _dats[not b,j+1].y); - t3:=abs(t2-_ergebnis[offset + i].x-_ergebnis[offset + i - 1].y+_ergebnis[offset + i - 1].x); + t2:=t2 + (_dats[not b,j+1]['x']-_dats[not b,j]['x']) * (_dats[not b,j]['y'] - t1) / (_dats[not b,j]['y'] - _dats[not b,j+1]['y']); + t3:=abs(t2-_ergebnis[offset + i]['x']-_ergebnis[offset + i - 1]['y']+_ergebnis[offset + i - 1]['x']); if (bestPos<0) or (bestErr>t3) then begin bestErr:=t3; bestPos:=t2; @@ -1000,19 +999,19 @@ begin inc(j); until j>_stops[not b]-_starts[not b]; end; - _ergebnis[offset + i].y:=bestPos; + _ergebnis[offset + i]['y']:=bestPos; end; end; for i:=_stops[false]-_starts[false]+1 to length(_ergebnis)-1 do begin - t1:=_ergebnis[i].x; - _ergebnis[i].x:=_ergebnis[i].y; - _ergebnis[i].y:=t1; + t1:=_ergebnis[i]['x']; + _ergebnis[i]['x']:=_ergebnis[i]['y']; + _ergebnis[i]['y']:=t1; end; for i:=0 to length(_ergebnis)-1 do begin - _ergebnis[i].y:=_ergebnis[i].y+_absShift; - t1:=(_ergebnis[i].y+_ergebnis[i].x)/2; - _ergebnis[i].y:=(_ergebnis[i].y-_ergebnis[i].x)/2; - _ergebnis[i].x:=t1; + _ergebnis[i]['y']:=_ergebnis[i]['y']+_absShift; + t1:=(_ergebnis[i]['y']+_ergebnis[i]['x'])/2; + _ergebnis[i]['y']:=(_ergebnis[i]['y']-_ergebnis[i]['x'])/2; + _ergebnis[i]['x']:=t1; end; fertig:=true; end; -- cgit v1.2.3-70-g09d2