program tkygras;                {Auswertung von Thermospannungsmessungen,}
                                {die mit TKMGRAF3 durchgef�hrt wurden}
uses  dos,crt,graph;            {numerische Ermittlung einer Funktion}
                                {alpha(T) als Stufenfunktion}
                                 {bisher beste L”sung}

const mxhi=120;myhi=120;hh1=15;hh2=15;myh=35;okmax=100;maxgrad=50;
var US,tem1,tem2,ur,dt,utm,ust:real;
    b1,b2,b4,b11,b12,b14,b35,b36:char;
    Name,name3:string;
    tex,f:Text;
    i,l,anz,imin,imax,imini,imaxi,b3,b5,b6,b7,b8,ii,b13,grad:Integer;
    z0,z1,z2,z3,z4,z5: real;

var graphdriver,graphmode:integer;
    mx,my,xp,yp,x0,y0,j,find,d,nr,find1,ok:integer;
    r,dx,dy,xl,yl,invdifx:real;
    xmin,xmax,x,y,ymin,ymax,difx,dify,invdify,deltax,deltay:real;
    sanz,strzahl:string;
    xr,yr,symin,symax,rr,xmini,xmaxi,temin,temax,dte:real;
    word1,word2,word3,word4,word5,sgrad,sok:string;
    c: array[1..maxgrad] of real;
    de: array[1..maxgrad] of extended;
    b:array[1..maxgrad,1..maxgrad] of extended;
    stelle:array [1..maxgrad] of integer;
    produkt:real;
    otem1,otem2,ous:array[-okmax..okmax] of real;
    oi,xmintxt,xmaxtxt,ymintxt,ymaxtxt,maxpix,nks,ischrift:integer;
    xmindef,xmaxdef,ymindef,ymaxdef:real;

procedure ig;
 begin
 initgraph(graphdriver,graphmode,'c:\sprachen\tp\bgi');
 end;

function kor(t_m1,t_m2:real):real;
 var ergebnis:real;
begin
if b6=1 then ergebnis:=0 else begin
 ergebnis:=z0*(t_m1-t_m2);
 ergebnis:=ergebnis+z1*(sqr(t_m1)-sqr(t_m2))/2;
 ergebnis:=ergebnis+z2*(t_m1*sqr(t_m1)-t_m2*sqr(t_m2))/3;
 ergebnis:=ergebnis+z3*(sqr(t_m1*t_m1)-sqr(t_m2*t_m2))/4;
 ergebnis:=ergebnis+z4*(t_m1*sqr(t_m1*t_m1)-t_m2*sqr(t_m2*t_m2))/5 ;
 ergebnis:=ergebnis+z5*(sqr(t_m1*t_m1*t_m1)-sqr(t_m2*t_m2*t_m2))/6 ;
end;
 kor:=ergebnis;
end;

function we(iwe:integer;twe:real):integer;
var wewe:integer;
begin
if twe<(temin+(iwe-1)*dte) then wewe:=0
else wewe:=1;
if iwe=1 then wewe:=1;
we:=wewe;
end;


function alpha100(xx:real):real;
var ialph,kalph:integer;
       yy:real;
begin
yy:=0;
for ialph:=1 to grad do begin
yy:=yy+c[ialph]*we(ialph,xx);end;
alpha100:=yy;
end;

function te(ite:integer):real;
begin
te:=temin+ite*dte;
end;

procedure matrixtrafo(j1:integer); {++++++++++++++++++++++++++++++++}
var bb:array[1..maxgrad] of real;
    dede:real;
    j2,jmax:integer;
begin
if (j1<grad) then begin
jmax:=j1;
for j2:=j1+1 to grad do begin
if (abs(b[jmax,j1])<abs(b[j2,j1])) then jmax:=j2;end;
for j2:=1 to grad do begin
bb[j2]:=b[j1,j2];b[j1,j2]:=b[jmax,j2];b[jmax,j2]:=bb[j2];end;
dede:=de[j1];de[j1]:=de[jmax];de[jmax]:=dede;
end;
end;

procedure koeff;
var detb:real;
    faktor:extended;
    ka,ki,i1,i2:integer;
begin
for i1:=1 to grad do begin
matrixtrafo(i1);
i2:=0;
repeat
inc(i2);
stelle[i1]:=i2;
until b[i1,i2]<>0;
 if i1<grad then begin
  for ka:=i1+1 to grad do begin
  faktor:=b[ka,i2]/b[i1,i2];
   for ki:=i2 to grad do begin
   b[ka,ki]:=b[ka,ki]-faktor*b[i1,ki];
   end;
  b[ka,i2]:=0;
  de[ka]:=de[ka]-faktor*de[i1];
  end;
 end;
end;
end;

procedure koeff2;
var i1,i2:integer;
begin
for i1:=grad downto 2 do begin
  c[stelle[i1]]:=de[i1]/b[i1,stelle[i1]];
    for i2:=i1-1 downto 1 do begin
    de[i2]:=de[i2]-c[stelle[i1]]*b[i2,stelle[i1]];
    end;
end;
c[stelle[1]]:=de[1]/b[1,stelle[1]];
end;

function komma(xxx:real):integer;
var il:integer;
begin
il:=0;
while xxx<1 do begin
xxx:=xxx*10;inc(il);end;
komma:=il;
end;

begin
repeat
 clrscr;
 writeln('Programm zur Auswertung von Thermospannungsmeáwerten, die mit');
 writeln('Hilfe des Programms tkmgraf3 aufgenommen wurden.');
 rr:=0;
 writeln('Eingabe des Namen der Datei, in der sich die Meáwerte befinden, die');
 write('ausgewertet werden sollen (mit Extension): ');
 readln(name);
 write('Probennummer: ');readln(word1);
 write('Schichtdicke in nm: ');readln(word2);
 write('Meádatum: ');readln(word3);
 write('Meástrom in mA: ');readln(word4);
 write('Target 1 (1) / Target 2 (2) ');readln(b7);
 write('Probe kristallin (1) / amorph (2) ');readln(b8);
{ writeln('Geben Sie die šberschrift der Grafik ein! (max. eine Zeile!)');
 readln(word5);}
 write('Auswertung relative Thermokraft (1) / absolute Thermokraft (2) ');
 readln(b6);
 write('Autografik (j/n) ');readln(b35);
   if b35='n' then begin
   xmindef:=0;xmaxdef:=300;ymindef:=-10;ymaxdef:=6;
{   write('Xmin : ');readln(xmindef);
   write('Xmax : ');readln(xmaxdef);
   write('Ymin : ');readln(ymindef);
   write('Ymax : ');readln(ymaxdef); }
   end;
if b35='j' then begin
 write('Soll der Koordinatenursprung in der Grafik enthalten sein? (j/n) ');
 readln (b4);end;
 write ('Gitter mitzeichnen? (j/n) ');readln(b36);
 utm:=0;ust:=0;
 writeln('Es wird alpha im betrachteten Temperaturbereich als Stufenfunktion');
 writeln('ermittelt. Wieviele Stufen sollen berechnet werden? (3<n<',maxgrad+1,')');
  repeat
 write('n = ');readln(grad);str(grad,sgrad);
 until ((3<grad) and (grad<=maxgrad));
 repeat
 write('Eingabe k-Wert (1...100): k = ');readln(ok);str(ok,sok);
 until ((0<ok)and(ok<101));
 write('Berechnung im gesamten Temperaturbereich ? (j/n) ');readln(b14);
 if b14<>'j' then begin
 writeln('Geben Sie den Temperaturbereich an, in dem die Stufenfunktionen berechnet');
 write('werden sollen!   Tmin/K:  ');readln(temin);
 write('                 Tmax/K:  ');readln(temax);
 end else begin
 write ('Nur Auswertung von Meáwerten mit T < 273 K ? (j/n) ');readln(b12);
 end;
 writeln('Eingabe des Namen der Datei, in der die Koeffizienten der Stufenfunktion');
 write('gespeichert werden: ');readln(name3);
 write('Sind alle Eingaben richtig? (j/n) ');readln(b11);
until b11='j';
 assign(f,'c:\daten\th_delta.txt');
 reset(f);
 readln(f,z0);
 readln(f,z1);
 readln(f,z2);
 readln(f,z3);
 readln(f,z4);
 readln(f,z5);
 close(f);
for i:=1 to maxgrad do begin
c[i]:=0;de[i]:=0;stelle[i]:=0;
for ii:=1 to maxgrad do b[i,ii]:=0;
end;

for b3:=2 to 4 do begin
 assign(tex,name);
 reset(tex);
if (b3<>3) then begin
repeat
readln(tex,anz);
readln(tex,tem1);
readln(tex,tem2);
readln(tex,us);us:=us-ust-kor(tem1,tem2);
readln(tex,dt);
readln(tex,rr);
xmin:=(tem1+tem2)/2;xmax:=xmin;
if b3=4 then begin ymin:=(tem1-tem2);ymax:=ymin;end;
until utm<=abs(us);
while (not eof(tex)) do begin
  readln(tex,anz);
  readln(tex,tem1);
  readln(tex,tem2);
  readln(tex,us);us:=us-ust-kor(tem1,tem2);
  readln(tex,dt);
  readln(tex,rr);
 if utm<=abs(us) then begin
  if (tem1+tem2)/2<xmin then xmin:=(tem1+tem2)/2;
  if (tem1+tem2)/2>xmax then xmax:=(tem1+tem2)/2;
 if b3=4 then begin
  if (tem1-tem2)<ymin then ymin:=(tem1-tem2);
  if (tem1-tem2)>ymax then ymax:=(tem1-tem2);
 end;end;
end;
end;
if b3<>3 then begin
if((b12='j')and(xmax>273)) then xmax:=273;
if ((b14='j')and(b3=2)) then begin temin:=xmin;temax:=xmax;end;
xmin:=xmin-0.1*abs(xmax-xmin);xmax:=xmax+0.09*abs(xmax-xmin);
end;
if b3<>2 then begin
ymin:=ymin-0.1*abs(ymax-ymin);ymax:=ymax+0.09*abs(ymax-ymin);
end;
if b4='j' then begin
xmin:=0;
if ymin>0 then ymin:=0;
if ymax<0 then ymax:=0;end;
if ((b3<>4)and(b35='n')) then begin
xmin:=xmindef;xmax:=xmaxdef;ymin:=ymindef;ymax:=ymaxdef;end;
if ((b3=4)and(b35='n')) then begin
xmin:=xmindef;xmax:=xmaxdef;end;


graphdriver:=detect;
ig;
mx:=getmaxx-20;my:=getmaxy-50;
    invdifx:=(mx-mxhi)/(xmax-xmin);
for i:=mxhi to mx do begin putpixel(i,myhi,1);
                           putpixel(i,my,1);end;
for i:=myhi to my do begin putpixel(mxhi,i,1);
                           putpixel(mx,i,1);end;

if b3<>2 then begin
dx:=(xmax-xmin);
if dx>0 then begin
r:=1;
    if dx>=100 then begin
           while dx>=100 do begin
                            dx:=dx/10;
                            r:=r*10;
                            end;                           {Ermittlung des}
                     end;                                  {X-Rasters}
    while dx<10 do begin
                   dx:=dx*10;
                   r:=r/10;
                   end;
    if dx>50 then deltax:=10*r
    else if dx>20 then deltax:=5*r
    else deltax:=2*r;
    nks:=komma(deltax);
    settextstyle(0,horizdir,0);
    settextjustify(centertext,righttext);
    invdifx:=(mx-mxhi)/(xmax-xmin);
    x:=int(xmin/deltax)*deltax;
    ischrift:=0;
    for j:=12 downto -1 do begin
    xr:=x+j*deltax;
    xl:=mxhi+(-xmin+x+j*deltax)*invdifx;
    x0:=round(xl);
    if (x0>=mxhi) and (x0<=mx) then begin
    inc(ischrift);
if ischrift<>2 then begin
str(xr:0:nks,strzahl);
outtextxy(x0,my+hh1,strzahl);
end else
outtextxy(x0,my+hh1,'T[K]');
if b36='j' then maxpix:=myhi else maxpix:=my-6;
    for i:=maxpix to my+6 do begin
    putpixel(x0,i,8);
    end;
    end;
    end;
x0:=round(mxhi-xmin*invdifx);
if (x0>=mxhi) and(x0<=mx) then begin
    for i:=myhi to my do begin
    putpixel(x0,i,4);
    end;
    end;
end;
{+++++++++++++++++++++++++++++}

xl:=mxhi+(-xmin+x+j*deltax)*invdifx;
dy:=ymax-ymin;
if dy>0 then begin
r:=1;
    if dy>=100 then begin
           while dy>=100 do begin
                            dy:=dy/10;
                            r:=r*10;
                            end;                           {Ermittlung des}
                     end;                                  {Y-Rasters}
    while dy<10 do begin
                   dy:=dy*10;
                   r:=r/10;
                   end;
    if dy>50 then deltay:=10*r
    else if dy>20 then deltay:=5*r
    else deltay:=2*r;
        nks:=komma(deltay);
        settextstyle(0,horizdir,0);                              
        settextjustify(righttext,centertext);
    y:=int(ymax/deltay)*deltay;
invdify:=(my-myhi)/(ymax-ymin);
    ischrift:=0;
    for j:=-1 to 12 do begin
    yr:=y-j*deltay;
    yl:=myhi+(ymax-y+j*deltay)*invdify;
    y0:=round(yl);
if((y0>=myhi)and(y0<=my)) then begin
    inc(ischrift);
    if ischrift<>2 then begin
str(yr:0:nks,strzahl);
outtextxy(mxhi-hh1,y0,strzahl);
end else begin
if b3=3 then outtextxy(mxhi-hh1,y0,'S[æV/K]');
if b3=4 then outtextxy(mxhi-hh1,y0,'T1-T2[K]');end;
if b36='j' then maxpix:=mx else maxpix:=mxhi+6;
    for i:=mxhi-6 to maxpix do begin                    
    putpixel(i,y0,8);
    end;
    end;
    end;

y0:=round(myhi+ymax*invdify);
if (y0>=myhi) and (y0<=my) then begin
    for i:=mxhi to mx do begin
    putpixel(i,y0,4);
    end; 
    end;
end;
end;

if b3=2 then begin
xmini:=(xmin+xmax)/2;xmaxi:=(xmin+xmax)/2;
end;
dte:=(temax-temin)/grad;
reset(tex);oi:=0;
while not eof(tex) do begin
    readln(tex,anz);
    readln(tex,tem1);
    readln(tex,tem2);
    readln(tex,us);us:=us-ust-kor(tem1,tem2);
    readln(tex,dt);
    readln(tex,rr);
for i:=-ok to ok-1 do begin
    otem1[i]:=otem1[i+1];
    otem2[i]:=otem2[i+1];
    ous[i]:=ous[i+1];end;
    otem1[ok]:=tem1;
    otem2[ok]:=tem2;
    ous[ok]:=us;
    inc(oi);
    x:=(tem1+tem2)/2;
if (((b12='j')and(tem1<273)and(tem2<273))or(b12<>'j'))then begin
    if utm<=abs(us) then begin
if b3=2 then begin
if tem1<xmini then xmini:=tem1;
if tem2<xmini then xmini:=tem2;
if tem1>xmaxi then xmaxi:=tem1;
if tem2>xmaxi then xmaxi:=tem2;
if ((tem1+tem2)/2>temin)and((tem1+tem2)/2<temax)and(oi>2*ok) then begin
tem1:=tem1;tem2:=tem2;us:=us;
for i:=1 to grad do begin
for ii:=1 to grad do begin
produkt:=((otem2[ok]-otem2[-ok])*we(i,otem2[0])-(otem1[ok]-otem1[-ok])*we(i,otem1[0]))*
((otem2[ok]-otem2[-ok])*we(ii,otem2[0])-(otem1[ok]-otem1[-ok])*we(ii,otem1[0]));
b[ii,i]:=b[ii,i]+produkt;
end;
produkt:=((otem2[ok]-otem2[-ok])*we(i,otem2[0])-(otem1[ok]-otem1[-ok])*we(i,otem1[0]));
produkt:=produkt*(ous[ok]-ous[-ok]);
de[i]:=de[i]+produkt;
end;
end;
end;
 if b3=4 then begin y:=(tem1-tem2);
                    x:=mxhi+(x-xmin)*invdifx;xp:=round(x);   { Meáwerte }
                    y:=myhi+(ymax-y)*invdify;yp:=round(y);
                    circle(xp,yp,2);
             end;end;
end;end;

if b3=2 then begin
   koeff;
   koeff2;
   sound(1000);delay(1000);nosound;
   if xmin>xmini then xmini:=xmin;
   if xmax<xmaxi then xmaxi:=xmax;
   if temin>xmini then xmini:=temin;
   if temax<xmaxi then xmaxi:=temax;
imin:=round(mxhi+(xmini-xmin)*invdifx);
imax:=round(mxhi+(xmaxi-xmin)*invdifx);
end;
if ((b3=2)or(b3=3)) then begin
   imini:=imin;
   imaxi:=imax;
   setcolor(13);
   for i:=imini to imaxi do begin
   x:=xmin+(i-mxhi)*(xmax-xmin)/(mx-mxhi);
   x:=x;
   y:=-1*alpha100(x);                           {Vorzeichenwechsel!}
if b3=2 then begin
if i=imini then begin ymin:=y;ymax:=y end else begin
if y<ymin then ymin:=y;if y>ymax then ymax:=y;end;
end;
if b3=3 then begin
   y:=myhi+(ymax-y)*invdify;yp:=round(y);
   if i>imini then lineto(i,yp);moveto(i,yp);
   putpixel(i,yp,13);
   end;
end;
   setcolor(15);
end;

if b3<>2 then begin
settextjustify(lefttext,centertext);
 settextstyle(0,horizdir,0);
outtextxy(15,75,'Probennr.: ');
outtextxy(110,75,word1);
outtextxy(round(mx/2)-20,75,'Schichtdicke/nm: ');
outtextxy(round(mx/2+130),75,word2);
outtextxy(15,45,'Meádatum: ');
outtextxy(110,45,word3);
outtextxy(round(mx/2)-20,45,'Meástrom/mA: ');
outtextxy(round(mx/2+130),45,word4);
if b6=1 then outtextxy(15,15,'Relative Thermokraft')
else outtextxy(15,15,'Absolute Thermokraft');
if b7=1 then outtextxy(265,15,'Target 1')
else outtextxy(265,15,'Target 2');
if b8=1 then outtextxy(415,15,'Probe kristallin')
else outtextxy(415,15,'Probe amorph');
{outtextxy(215,15,word5);}
if b3=3 then begin
outtextxy(round(mx*3/4)+120,45,sgrad);
outtextxy(round(mx*3/4)+120,75,sok);
settextjustify(righttext,centertext);
outtextxy(round(mx*3/4)+120,45,'n = ');
outtextxy(round(mx*3/4)+120,75,'k = ');
end;
end;

close(tex);
if b3<>2 then readln;
closegraph;
end;
writeln;
 assign(tex,name3);
 rewrite(tex);
for i:=1 to grad do begin
x:=temin+(i-0.5)*dte;
y:=-alpha100(x);
writeln(tex,x);
writeln(tex,y);
end;
close(tex);
writeln('Programmende.');
readln;
end.