macheKomplex und istKomplex neu

1 files changed, 164 insertions, 39 deletions

diff --git a/epostunit.pas b/epostunit.pas
index ed65ab1..c393734 100644
--- a/epostunit.pas
+++ b/epostunit.pas
@@ -5,7 +5,8 @@ unit epostunit;
 interface
 
 uses
-  classes, sysutils, mystringlistunit, werteunit, typenunit, process, lowlevelunit, matheunit, fftunit;
+  classes, sysutils, mystringlistunit, werteunit, typenunit, process, lowlevelunit,
+  matheunit, fftunit, randomunit;
 
 type
   tBmpHeader = packed record
@@ -117,6 +118,7 @@ type
     procedure wMinW(miW: extended);
     function rMaxW: extended;
     procedure wMaxW(maW: extended);
+    function istKomplex: boolean;
     function xScale: extended;
     function tScale: extended;
     function dichtenParameterErkannt(sT: boolean; s: string; var bekannteBefehle: tMyStringList; threads,xMin,xMax,tMin,tMax: longint): boolean;
@@ -159,6 +161,7 @@ type
     function entferneArtefakte(sT: boolean; f: tMyStringList; threads: longint): boolean;
     function extrahiereEinhuellende(sT: boolean; f: tMyStringList; threads: longint; Warn: tWarnstufe): boolean;
     function extrahierePhase(sT: boolean; f: tMyStringList; threads: longint; Warn: tWarnstufe): boolean;
+    function macheKomplex(sT: boolean; f: tMyStringList; threads: longint): boolean;
     function berechneIntegral(sT: boolean; f: tMyStringList; threads: longint; quelle: tWerte): boolean;
     function berechneFFT(sT: boolean; f: tMyStringList; threads: longint; Warn: tWarnstufe): boolean;
     function berechneFFT2d(sT: boolean; f: tMyStringList; threads: longint; Warn: tWarnstufe): boolean;
@@ -169,7 +172,7 @@ type
     procedure spiegle(threads: longint); overload;
     procedure spiegle(threads,tMin,tMax: longint); overload;
     procedure fuelleMitDummys(sT: boolean);
-    procedure verschiebe(threads: longint; richtung: tIntPoint; komplex: boolean);
+    procedure verschiebe(threads: longint; richtung: tIntPoint);
     procedure ermittlePhasenWinkel(threads: longint);
     procedure fft2dNachbearbeitung(threads: longint; nB: tFFTDatenordnung);
     procedure schreibeWert(var f: textfile; x,y: longint);
@@ -245,16 +248,17 @@ type
     nummer,mT,breite,wBreite,hoehe,
       wHoehe,gesBreite,lOf,oOf,rOf,uOf: longint;
     ws:                                 tWerteArray;
-    xMi,xMa,tMi,tMa,xPMi,xPMa,komplex:  longint;
+    xMi,xMa,tMi,tMa,xPMi,xPMa:          longint;
     xZ,yZ:                              extended;
     nbs:                                tTransformationArray;
     farben:                             tRGBArray;
+    istKomplex:                         array of boolean;
     wertes:                             array of tExtendedArray;
     anzahlens:                          array of tLongintArray;
     pals:                               array of tPalette;
     rahmen:                             boolean;
     beschr:                             pTBeschriftungen;
-    constructor create(i,maxThreads,iBreite,iHoehe,lO,oO,rO,uO: longint; const wes: tWerteArray; xMin,xMax,tMin,tMax,kmpl: longint; xZoom,yZoom: extended; nachbearbeitungen: tTransformationArray; paletten: tPalettenArray; beschri: pTBeschriftungen; rm: boolean);
+    constructor create(i,maxThreads,iBreite,iHoehe,lO,oO,rO,uO: longint; const wes: tWerteArray; xMin,xMax,tMin,tMax: longint; xZoom,yZoom: extended; nachbearbeitungen: tTransformationArray; paletten: tPalettenArray; beschri: pTBeschriftungen; rm: boolean);
     destructor destroy; override;
     procedure stExecute; override;
     procedure initAnzahlensFuerKontur;
@@ -389,8 +393,7 @@ type
     we:              tWerte;
     xMi,xMa,tMi,tMa: longint;
     rtg:             tIntPoint;
-    kmpl:            boolean;
-    constructor create(werte: tWerte; xMin,xMax,tMin,tMax: longint; richtung: tIntPoint; komplex: boolean);
+    constructor create(werte: tWerte; xMin,xMax,tMin,tMax: longint; richtung: tIntPoint);
     procedure stExecute; override;
   end;
   tPhasenWinkelThread = class(tLogThread)
@@ -399,6 +402,14 @@ type
     constructor create(werte: tWerte; xMin,xMax: longint);
     procedure stExecute; override;
   end;
+  tKomplexMachThread = class(tLogThread)
+    we:      tWerte;
+    xMi,xMa: longint;
+    kmm:     tKomplexMachModus;
+    mT:      tMersenneTwister;
+    constructor create(werte: tWerte; xMin,xMax: longint; komplexMachModus: tKomplexMachModus; randomSeed: longword);
+    procedure stExecute; override;
+  end;
 
 function neuePalette(f: tMyStringList): boolean;
 function initBmpHeader(w,h: longint): tBmpHeader;
@@ -933,6 +944,18 @@ begin
   end{of case};
 end;
 
+function tWerte.istKomplex: boolean;
+begin
+  case genauigkeit of
+    gSingle:
+      result:=sWerte.params.istKomplex;
+    gDouble:
+      result:=dWerte.params.istKomplex;
+    gExtended:
+      result:=eWerte.params.istKomplex;
+  end{of case};
+end;
+
 function tWerte.findeAlleDateien(nam: string; var dat: tGenerischeInputDateiInfoArray; vorlage: tGenerischeInputDateiInfo): boolean;
 var
   err: longint;
@@ -3921,7 +3944,7 @@ begin
   betraege.free;
 
   gibAus('zentrieren '+tIntPointToStr(maxPos),3);
-  verschiebe(threads,maxPos,true);
+  verschiebe(threads,maxPos);
 
   gibAus('berechne inverse x-FFT ...',3);
   if not fft(threads,false,true,doAlleResIms,doAlleResIms,nil,0,pvFehler,wsLasch) then begin
@@ -3943,6 +3966,80 @@ begin
   result:=true;
 end;
 
+function tWerte.macheKomplex(sT: boolean; f: tMyStringList; threads: longint): boolean;
+var
+  bekannteBefehle:  tMyStringList;
+  komplexMachModus: tKomplexMachModus;
+  kmThreads:        array of tKomplexMachThread;
+  s:                string;
+  i:                longint;
+  fertig:           boolean;
+begin
+  result:=false;
+  if istKomplex then begin
+    gibAus('Werte '''+bezeichner+''' sind bereits komplex!',3);
+    exit;
+  end;
+
+  komplexMachModus:=kmmImNull;
+
+  bekannteBefehle:=tMyStringList.create;
+  repeat
+    if not f.metaReadln(s,true) then begin
+      gibAus('Unerwartetes Dateiende!',3);
+      exit;
+    end;
+    bekannteBefehle.clear;
+    if istDasBefehl('Ende',s,bekannteBefehle,false) then break;
+    if istDasBefehl('rein reell',s,bekannteBefehle,false) then begin
+      komplexMachModus:=kmmImNull;
+      continue;
+    end;
+    if istDasBefehl('rein imaginär',s,bekannteBefehle,false) then begin
+      komplexMachModus:=kmmReNull;
+      continue;
+    end;
+    if istDasBefehl('zufällige Phase',s,bekannteBefehle,false) then begin
+      komplexMachModus:=kmmPhZuf;
+      continue;
+    end;
+    bekannteBefehle.sort;
+    gibAus('Verstehe Option '''+s+''' nicht beim Komplexmachen!'#10'Ich kenne:'#10+bekannteBefehle.text,3);
+    bekannteBefehle.free;
+    exit;
+  until false;
+
+  _tSiz:=_tSiz*2;
+  if sT then begin
+    result:=true;
+    exit;
+  end;
+
+  holeRAM(3);
+
+  setLength(kmThreads,threads);
+  for i:=0 to length(kmThreads)-1 do
+    kmThreads[i]:=tKomplexMachThread.create(
+      self,
+      round(i/threads*(_tSiz div 2)),
+      round((i+1)/threads*(_tSiz div 2))-1,
+      komplexMachModus,
+      random(high(longword))
+    );
+
+  repeat
+    fertig:=true;
+    for i:=0 to length(kmThreads)-1 do
+      fertig:=fertig and kmThreads[i].fertig;
+    if not fertig then sleep(10);
+  until fertig;
+
+  for i:=0 to length(kmThreads)-1 do
+    kmThreads[i].free;
+
+  result:=true;
+end;
+
 function tWerte.berechneIntegral(sT: boolean; f: tMyStringList; threads: longint; quelle: tWerte): boolean;
 var
   i,tMin,tMax,xMin,xMax: longint;
@@ -4238,7 +4335,7 @@ var
   nachbearbeitungen:                     tTransformationArray;
   Ausschnitt:                            tTransformation;
   bilderThreads:                         array of tBilderThread;
-  fertig,rahmen,komplex:                 boolean;
+  fertig,rahmen:                         boolean;
   img:                                   file;
   achsen:                                array of tAchse;
   fontRenderer:                          tFontRenderer;
@@ -4302,7 +4399,6 @@ begin
   Ausschnitt:=nil;
   rahmen:=false;
   fontRenderer:=nil;
-  komplex:=false;
   bekannteBefehle:=tMyStringList.create;
   repeat
     if not f.metaReadln(s,true) then begin
@@ -4369,10 +4465,6 @@ begin
       rahmen:=true;
       continue;
     end;
-    if istDasBefehl('komplex',s,bekannteBefehle,false) then begin
-      komplex:=true;
-      continue;
-    end;
     if istDasBefehl('Nachbearbeitung:',s,bekannteBefehle,true) then begin
       if s[1] in ['0'..'9'] then
         i:=strtoint(erstesArgument(s))
@@ -4393,7 +4485,7 @@ begin
       continue;
     end;
     if istDasBefehl('Legende:',s,bekannteBefehle,true) then begin
-      if komplex then begin
+      if istKomplex then begin
         gibAus('Ich kann keine Legende zu komplexen Werten erstellen!',3);
         aufraeumen;
         exit;
@@ -4519,7 +4611,7 @@ begin
   xMin:=max(xMin,0);
   xMax:=min(xMax,_xSteps-1);
   tMin:=max(tMin,0);
-  tMax:=min(tMax,_tSiz div (1+byte(komplex))-1);
+  tMax:=min(tMax,_tSiz div (1+byte(istKomplex))-1);
   gibAus('('+intToStr(xMin)+'-'+intToStr(xMax)+'x'+intToStr(tMin)+'-'+intToStr(tMax)+')',3);
   gibAus(' ('+floatToStr(transformationen.xStart)+'-'+floatToStr(transformationen.xStop)+' x '+floatToStr(transformationen.tStart)+'-'+floatToStr(transformationen.tStop)+')',3);
   breite:=round((xMax-xMin+1)*xZoom);
@@ -4595,7 +4687,7 @@ begin
 
   setLength(bilderThreads,maxThreads);
   for i:=0 to length(bilderThreads)-1 do
-    bilderThreads[i]:=tBilderThread.create(i,length(bilderThreads),breite,hoehe,lOf,oOf,rOf,uOf,quellen,xMin,xMax,tMin,tMax,byte(komplex)*_tSiz div 2,xZoom,yZoom,nachbearbeitungen,paletten,@beschriftungen,rahmen);
+    bilderThreads[i]:=tBilderThread.create(i,length(bilderThreads),breite,hoehe,lOf,oOf,rOf,uOf,quellen,xMin,xMax,tMin,tMax,xZoom,yZoom,nachbearbeitungen,paletten,@beschriftungen,rahmen);
   for i:=0 to length(bilderThreads)-1 do begin
     gibAus('starte Thread '+intToStr(i)+' ...',1);
     bilderThreads[i].suspended:=false;
@@ -5014,7 +5106,7 @@ begin
     end;
 end;
 
-procedure tWerte.verschiebe(threads: longint; richtung: tIntPoint; komplex: boolean);
+procedure tWerte.verschiebe(threads: longint; richtung: tIntPoint);
 var
   einheitsZelle:     tIntPoint;
   teilRichtung:      char;
@@ -5022,7 +5114,7 @@ var
   i:                 longint;
   fertig:            boolean;
 begin
-  einheitsZelle:=berechneEinheitsZelle(richtung,intPoint(_xSteps,_tSiz div (1+byte(komplex))));
+  einheitsZelle:=berechneEinheitsZelle(richtung,intPoint(_xSteps,_tSiz div (1+byte(istKomplex))));
 
   teilRichtung:=char(ord('x')+byte(einheitsZelle['y']>einheitsZelle['x']));
   if einheitsZelle[teilRichtung]<threads then
@@ -5038,8 +5130,7 @@ begin
           round((i+1)/threads*einheitsZelle['x']-1),
           0,
           einheitsZelle['y']-1,
-          richtung,
-          komplex
+          richtung
         )
     else
       verschiebeThreads[i]:=
@@ -5049,8 +5140,7 @@ begin
           einheitsZelle['x']-1,
           round(i/threads*einheitsZelle['y']),
           round((i+1)/threads*einheitsZelle['y']-1),
-          richtung,
-          komplex
+          richtung
         );
   repeat
     sleep(10);
@@ -5681,9 +5771,9 @@ end;
 
 // tBilderThread ***************************************************************
 
-constructor tBilderThread.create(i,maxThreads,iBreite,iHoehe,lO,oO,rO,uO: longint; const wes: tWerteArray; xMin,xMax,tMin,tMax,kmpl: longint; xZoom,yZoom: extended; nachbearbeitungen: tTransformationArray; paletten: tPalettenArray; beschri: pTBeschriftungen; rm: boolean);
+constructor tBilderThread.create(i,maxThreads,iBreite,iHoehe,lO,oO,rO,uO: longint; const wes: tWerteArray; xMin,xMax,tMin,tMax: longint; xZoom,yZoom: extended; nachbearbeitungen: tTransformationArray; paletten: tPalettenArray; beschri: pTBeschriftungen; rm: boolean);
 var
-  ii: longint;
+  ii,anzKompl: longint;
 begin
   inherited create;
   beschr:=beschri;
@@ -5707,7 +5797,6 @@ begin
   for ii:=0 to length(pals)-1 do
     pals[ii]:=paletten[ii];
   rahmen:=rm;
-  komplex:=kmpl;
   xPMi:=nummer*((gesBreite+lOf+rOf) div mT)-lOf;
   if nummer=mT-1 then xPMa:=gesBreite+rOf-1                 // der letzte Thread bekommt den Rest
   else xPMa:=(nummer+1)*((gesBreite+lOf+rOf) div mT)-lOf-1;
@@ -5717,14 +5806,19 @@ begin
 //  breite:=wBreite+lOf*byte(nummer=0)+rOf*byte(nummer=mT-1);
   hoehe:=oOf+wHoehe+uOf;
   gibAus('Werte: '+intToStr(xMi)+'-'+intToStr(xMa)+'x'+intToStr(tMi)+'-'+intToStr(tMa)+' ('+intToStr(wBreite)+'), Pixel: '+intToStr(xPMi)+'-'+intToStr(xPMa)+' ('+intToStr(breite)+')',1);
+  anzKompl:=0;
+  for ii:=0 to length(wes)-1 do
+    anzKompl:=anzKompl+1+byte(wes[ii].istKomplex);
   gibAus('Thread '+intToStr(nummer)+': hole '
-    +intToStr(round(((sizeOf(extended)+sizeOf(longint))*length(wes)*(1+byte(komplex))*wHoehe*wBreite+
+    +intToStr(round(((sizeOf(extended)+sizeOf(longint))*anzKompl*wHoehe*wBreite+
                       sizeOf(tRGB)*hoehe*breite)/1024/1024))+'MB RAM ...',1);
+  setLength(istKomplex,length(wes));
   setLength(wertes,length(wes));
   setLength(anzahlens,length(wes));
   for ii:=0 to length(wertes)-1 do begin
+    istKomplex[ii]:=wes[ii].istKomplex;
     setLength(anzahlens[ii],wHoehe*wBreite);
-    setLength(wertes[ii],(1+byte(komplex>0))*length(anzahlens[ii]));
+    setLength(wertes[ii],(1+byte(wes[ii].istKomplex))*length(anzahlens[ii]));
   end;
   setLength(farben,hoehe*breite);
   gibAus('Thread '+intToStr(nummer)+' hat jetzt seinen Speicher.',1);
@@ -5747,6 +5841,7 @@ end;
 procedure tBilderThread.stExecute;
 var
   i,j,k:   longint;
+  hLen:    int64;
   wert:    extended;
   b:       boolean;
   frb:     tRGB;
@@ -5754,9 +5849,9 @@ begin
   gibAus('Thread '+intToStr(nummer)+' gestartet!',1);
   for i:=0 to length(ws)-1 do
     case ws[i].genauigkeit of
-      gSingle:   if not ws[i].sWerte.zuPixelWerten(wHoehe,wBreite,xPMi,xMi,tMi,komplex,xZ,yZ,@wertes[i],@anzahlens[i]) then exit;
-      gDouble:   if not ws[i].dWerte.zuPixelWerten(wHoehe,wBreite,xPMi,xMi,tMi,komplex,xZ,yZ,@wertes[i],@anzahlens[i]) then exit;
-      gExtended: if not ws[i].eWerte.zuPixelWerten(wHoehe,wBreite,xPMi,xMi,tMi,komplex,xZ,yZ,@wertes[i],@anzahlens[i]) then exit;
+      gSingle:   if not ws[i].sWerte.zuPixelWerten(wHoehe,wBreite,xPMi,xMi,tMi,xZ,yZ,@wertes[i],@anzahlens[i]) then exit;
+      gDouble:   if not ws[i].dWerte.zuPixelWerten(wHoehe,wBreite,xPMi,xMi,tMi,xZ,yZ,@wertes[i],@anzahlens[i]) then exit;
+      gExtended: if not ws[i].eWerte.zuPixelWerten(wHoehe,wBreite,xPMi,xMi,tMi,xZ,yZ,@wertes[i],@anzahlens[i]) then exit;
     end{of case};
   for i:=0 to length(wertes)-1 do
     for j:=0 to length(anzahlens[i])-1 do
@@ -5771,9 +5866,10 @@ begin
         farben[i-xPMi+(j+oOf)*breite].rgbRed:=0;
         farben[i-xPMi+(j+oOf)*breite].rgbGreen:=0;
         farben[i-xPMi+(j+oOf)*breite].rgbBlue:=0;
-        for k:=0 to length(wertes)-1 do
-          for b:=false to komplex>0 do begin
-            wert:=(wertes[k,i-max(0,xPMi)+(j+byte(b)*komplex)*wBreite]/anzahlens[k,i-max(0,xPMi)+j*wBreite]-ws[k]._minW)/(ws[k]._maxW-ws[k]._minW);
+        for k:=0 to length(wertes)-1 do begin
+          hLen:=length(wertes[k]) div 2;
+          for b:=false to istKomplex[k] do begin
+            wert:=(wertes[k,i-max(0,xPMi)+(j+byte(b and istKomplex[k])*hLen)*wBreite]/anzahlens[k,i-max(0,xPMi)+j*wBreite]-ws[k]._minW)/(ws[k]._maxW-ws[k]._minW);
             wert:=max(0,min(wert,1));
             wert:=nbs[k].transformiereWert(wert);
             frb:=pals[k].wertZuFarbe(wert,b);
@@ -5783,6 +5879,7 @@ begin
               rgbBlue:= rgbBlue  + frb.rgbBlue;
             end;
           end;
+        end;
       end
       else begin
         farben[i-xPMi+(j+oOf)*breite].rgbRed:=
@@ -7147,7 +7244,7 @@ end;
 
 // tVerschiebeThread ***********************************************************
 
-constructor tVerschiebeThread.create(werte: tWerte; xMin,xMax,tMin,tMax: longint; richtung: tIntPoint; komplex: boolean);
+constructor tVerschiebeThread.create(werte: tWerte; xMin,xMax,tMin,tMax: longint; richtung: tIntPoint);
 begin
   inherited create;
   we:=werte;
@@ -7156,8 +7253,7 @@ begin
   tMi:=tMin;
   tMa:=tMax;
   rtg:=richtung;
-  kmpl:=komplex;
-  gibAus('VerschiebeThread erzeugt: '+intToStr(xMin)+'-'+intToStr(xMax)+' '+intToStr(tMin)+'-'+intToStr(tMax)+' '+tIntPointToStr(richtung)+' '+intToStr(byte(komplex)),1);
+  gibAus('VerschiebeThread erzeugt: '+intToStr(xMin)+'-'+intToStr(xMax)+' '+intToStr(tMin)+'-'+intToStr(tMax)+' '+tIntPointToStr(richtung)+' '+intToStr(byte(we.istKomplex)),1);
   suspended:=false;
 end;
 
@@ -7166,11 +7262,11 @@ begin
   gibAus('VerschiebeThread gestartet',1);
   case we.genauigkeit of
     gSingle:
-      we.sWerte.verschiebe(rtg,xMi,xMa,tMi,tMa,kmpl);
+      we.sWerte.verschiebe(rtg,xMi,xMa,tMi,tMa);
     gDouble:
-      we.dWerte.verschiebe(rtg,xMi,xMa,tMi,tMa,kmpl);
+      we.dWerte.verschiebe(rtg,xMi,xMa,tMi,tMa);
     gExtended:
-      we.eWerte.verschiebe(rtg,xMi,xMa,tMi,tMa,kmpl);
+      we.eWerte.verschiebe(rtg,xMi,xMa,tMi,tMa);
   end{of case};
   gibAus('VerschiebeThread beendet',1);
 end;
@@ -7201,6 +7297,35 @@ begin
   gibAus('PhasenWinkelThread beendet',1);
 end;
 
+// tKomplexMachThread **********************************************************
+
+constructor tKomplexMachThread.create(werte: tWerte; xMin,xMax: longint; komplexMachModus: tKomplexMachModus; randomSeed: longword);
+begin
+  inherited create;
+  we:=werte;
+  xMi:=xMin;
+  xMa:=xMax;
+  kmm:=komplexMachModus;
+  mT:=tMersenneTwister.create;
+  mT.init(randomSeed);
+  gibAus('KomplexMachThread erzeugt: '+intToStr(xMin)+'-'+intToStr(xMax),1);
+  suspended:=false;
+end;
+
+procedure tKomplexMachThread.stExecute;
+begin
+  gibAus('KomplexMachThread gestartet',1);
+  case we.genauigkeit of
+    gSingle:
+      we.sWerte.macheKomplex(xMi,xMa,kmm,mT);
+    gDouble:
+      we.dWerte.macheKomplex(xMi,xMa,kmm,mT);
+    gExtended:
+      we.eWerte.macheKomplex(xMi,xMa,kmm,mT);
+  end{of case};
+  gibAus('KomplexMachThread beendet',1);
+end;
+
 // sonstiges *******************************************************************
 
 function findePalette(out palette: tPalette; name: string): boolean;