<*+M2EXTENSIONS *>
<*HEAPLIMIT="1500000000"*>
<*-CHECKDIV *>
<*-CHECKRANGE *>
<*-GENFRAME*>
<*-COVERFLOW *>
<*-IOVERFLOW*>
<*-GENHISTORY*>
<*-GENDEBUG*>
<*-LINENO*>
<*-CHECKINDEX*>
<*-CHECKDINDEX*>
<*-NOPTRALIAS*>
<*-DOREORDER*>
<*-PROCINLINE*>
<*-GENPTRINIT*>
<*+STORAGE *>
<*CPU="PENTIUM"*>
<* IF __GEN_C__ THEN *>
  <*+M2EXTENSIONS *>
  <*+STORAGE *>
  <*-GENCTYPES*>
  <*+COMMENT*>
  <*-GENHISTORY*>
  <*-GENDEBUG*>
  <*-GENDATE*>
  <*-LINENO*>
  <*-CHECKINDEX*>
  <*-CHECKDNDEX*>
  <*+GENCDIV*>
  <*-GENKRC*>
  <*+NOOPTIMIZE*>
  <*-GENSIZE*>
  <*-ASSERT*>
  <*-CHECKNIL*>
  <*-COVERFLOW*>
  <*-IOVERFLOW*>
  <*-CHECKRANGE*>
  <*-CHECKSET*>
  <*-CHECKDIV*>
  <*-GENCONSTENUM*>
<* END *>

MODULE profile; (* altitude profile for radio link with srtm data by oe5dxl *)

FROM SYSTEM IMPORT FILL, ADR, INT16, INT8, CAST, MOVE, CARD8, CARD16, SHIFT;

FROM osi IMPORT Werr, WrStr, WrStrLn, WrInt, File, OpenWrite, OpenRead, Close, WrBin, 
     WrFixed, ALLOCATE, DEALLOCATE, NextArg, RdBin;

FROM math IMPORT sqrt, sin, cos, atan, exp, log;

FROM libsrtm IMPORT getsrtmlong, srtmdir, closesrtmfile, ATTRWATER, ATTRWOOD, ATTRURBAN;

FROM aprspos IMPORT distance, RAD, azimuth, EARTH (*, wgs84r, wgs84s*);
FROM aprsstr IMPORT POSITION, IntToStr, Append, Assign, DateToStr, FixToStr, Length,
     loctopos, TimeToStr, StrToTime, TIME, cleanfilename, InStr, StrToFix, StrCmp, StrToInt;

IMPORT pngwrite;

FROM imagetext IMPORT pIMAGE, writestr, strsize, fontsizex, fontsizey;

CONST
      FRAMEXR=5;
      MINM=20;
      NOALT=20000;
      PI=3.1415926535;
      PI2=PI*2.0;
      LF=12C;
      ERRORFONT=10;
      MAXAIRSHADOW=0.5;                   (* full red dB *)
      MINAIRSHADOW=0.1;                   (* begin red dB *)
      MINDIST=0.2;
      DEFMAXDIST=600.0;
      LEFTSPACE=3;                        (* begin of texts from left image margin *)
      FRESPERC=0.6;                       (* kernel diameter of fresnel zone *)

TYPE
     SET8=SET OF [0..7];
     SET16=SET OF [0..15];
     SET32=SET OF [0..31];
     COLS=(cEARTH, cHEAVENup, cHEAVENdown, cFRESIN, cFRESOUT, cOPTALT, cSCALERS, cMLINES, cFRESLINE,
           cMEADOW, cTEXT1, cTEXT2, cTEXT3, cTEXTSCALE, cBRANCH, cTRUNC, cWATER, cWOOD, cURBAN);

     POSITIONL=RECORD
                long,
                lat:LONGREAL;
              END;

(*
     PIX=RECORD r,g,b:CARD16 END;
     IMAGELINE=ARRAY OF PIX;
     IMAGE=ARRAY OF IMAGELINE;
     pIMAGE=POINTER TO IMAGE;
*)
     PATH=RECORD
       pos:POSITIONL;
       optalt,
       fresm,
       refrm, 
       zero,
       wood, 
       alt:LONGREAL;
       resol:REAL;
       attr:CARD8;
       treepri:BOOLEAN; 
     END;
     pPATH=POINTER TO ARRAY OF PATH;

VAR
    image     :pIMAGE;
    csvfn,
    imagefn   :ARRAY[0..1023] OF CHAR;
    xsize, ysize, linksize, fontx, fonty, fonttyp, framexl, frameyu, frameyd:INTEGER; 
    posa, posb:POSITIONL;
    alta, altb, anta, antb, refraction, mhz, dist, igamma, treesize, maxdist, kernel:LONGREAL;
    path      : pPATH;
    gammatab  : ARRAY[0..1023] OF CHAR;
    labela,
    labelb    : ARRAY[0..99] OF CHAR;
    realwood,
    opt,
    treedrawn : BOOLEAN;
    colours   : ARRAY[MIN(COLS)..MAX(COLS)] OF RECORD r,g,b:CARDINAL END;

PROCEDURE Error(text:ARRAY OF CHAR);
BEGIN
  Werr(text); Werr(" error abort"+LF);
  HALT
END Error;

PROCEDURE sqr(x:LONGREAL):LONGREAL; BEGIN RETURN x*x END sqr;

PROCEDURE posinval(VAR pos:POSITIONL);
BEGIN pos.long:=0.0; pos.lat:=0.0 END posinval;

PROCEDURE posvalid(pos-:POSITIONL):BOOLEAN;
BEGIN RETURN (pos.lat<>0.0) OR (pos.long<>0.0) END posvalid;


PROCEDURE wgs84s(lat, long, nn:LONGREAL; VAR x,y,z:LONGREAL);   (* km *)
VAR h, c:LONGREAL;
BEGIN
  h:=nn+EARTH;
  z:=h*sin(lat);
  c:=cos(lat);
  y:=h*sin(long)*c;
  x:=h*cos(long)*c;
END wgs84s;

PROCEDURE wgs84r(x,y,z:LONGREAL; VAR lat, long, heig:LONGREAL);   (* km *)
VAR h:LONGREAL;
BEGIN
  h:=x*x + y*y;
  IF ABS(x)>ABS(y) THEN
    long:=atan(y/x);
    IF x<0.0 THEN
      IF y>0.0 THEN long:=PI+long ELSE long:=long-PI END;
    END;
  ELSE
    long:=PI*0.5-atan(x/y);
    IF y<0.0 THEN long:=long-PI END;
  END;
  lat:=atan(z/sqrt(h));
  heig:=sqrt(h + z*z) - EARTH;
END wgs84r;


PROCEDURE StrToFixl(VAR x:LONGREAL; s-:ARRAY OF CHAR):BOOLEAN;
VAR r:REAL;
BEGIN
  IF NOT StrToFix(r, s) THEN RETURN FALSE END;
  x:=VAL(LONGREAL, r);
  RETURN TRUE
END StrToFixl; 

PROCEDURE makegammatab;
VAR c:CARDINAL;
    g:LONGREAL;
BEGIN
  g:=1.0/igamma;
  gammatab[0]:=0C;
  FOR c:=1 TO HIGH(gammatab) DO gammatab[c]:=CHR(TRUNC(exp(log(FLOAT(c)/1024.0)*g)*255.5)) END;
END makegammatab;


PROCEDURE wrpng;

  PROCEDURE pngc(c:CARD16):CARD8;
  BEGIN IF c<=HIGH(gammatab) THEN RETURN ORD(gammatab[c]) ELSE RETURN 255 END END pngc;

VAR pngimg:pngwrite.PNGPIXMAP;
    x, y, ret:INTEGER;
BEGIN
  ALLOCATE(pngimg.image, xsize*ysize*3);
  IF pngimg.image<>NIL THEN
    makegammatab;
    FOR y:=0 TO ysize-1 DO
      FOR x:=0 TO xsize-1 DO
        WITH image^[x][ysize-1-y] DO
          WITH pngimg.image^[x+y*xsize] DO red:=pngc(r); green:=pngc(g); blue:=pngc(b) END;
        END;
      END;
    END;
    pngimg.width:=xsize;
    pngimg.height:=ysize;
    ret:=pngwrite.writepng(imagefn, pngimg);
    DEALLOCATE(pngimg.image, xsize*ysize*3);
  ELSE Werr("png write out of memory"+LF) END;
END wrpng;


PROCEDURE background;            (* 20:80:220 *)
VAR i, y, rr, gg, bb:INTEGER;
BEGIN
  FOR y:=0 TO ysize-1 DO 
    rr:=(VAL(INTEGER, colours[cHEAVENup].r)*y + 
         VAL(INTEGER, colours[cHEAVENdown].r)*(ysize-y)) DIV ysize;
    gg:=(VAL(INTEGER, colours[cHEAVENup].g)*y +
         VAL(INTEGER, colours[cHEAVENdown].g)*(ysize-y)) DIV ysize;
    bb:=(VAL(INTEGER, colours[cHEAVENup].b)*y +
         VAL(INTEGER, colours[cHEAVENdown].b)*(ysize-y)) DIV ysize;
    FOR i:=0 TO xsize-1 DO
      WITH image^[i][y] DO
        r:=rr;
        g:=gg;
        b:=bb;
      END;
    END;
  END;
END background;

PROCEDURE errorimg(errmsg-:ARRAY OF CHAR);
VAR x:INTEGER;
BEGIN
  IF image=NIL THEN
    NEW(image, xsize, ysize);
    IF image=NIL THEN Error("out of memory") END;

  END;
  x:=(xsize-VAL(INTEGER,strsize(ERRORFONT, errmsg))) DIV 2;
  IF x<0 THEN x:=0 END;
  background;
  writestr(image, x, (ysize-VAL(INTEGER,fontsizey(ERRORFONT))) DIV 2, ERRORFONT, 0,0, 900,100,0, errmsg);
  wrpng;
  Error(errmsg);
END errorimg;

PROCEDURE card(s-:ARRAY OF CHAR; VAR p, res:CARDINAL; len:CARDINAL; VAR ok:BOOLEAN);
VAR n:CARDINAL;
BEGIN
  ok:=FALSE;
  n:=0;
  WHILE (p<len) & (s[p]=" ") DO INC(p) END;
  WHILE (p<len) & (s[p]>="0") & (s[p]<="9") DO
    n:=n*10+ORD(s[p])-ORD("0");
    ok:=TRUE;
    INC(p);
  END;
  IF n>1023 THEN ok:=FALSE END;
  IF ok THEN res:=n END;
END card;

PROCEDURE readcolours(fn:ARRAY OF CHAR);
VAR fd, len:INTEGER;
    p, lc:CARDINAL;
    ci:COLS;
    ok:BOOLEAN;
    s:ARRAY[0..5000] OF CHAR;
    h:ARRAY[0..100] OF CHAR;
BEGIN

  fd:=OpenRead(fn);
  IF fd<0 THEN s:="["; Append(s, fn); Append(s, "] colour file not found"); Error(s); END;
  len:=RdBin(fd, s, SIZE(s)-1);
  Close(fd);
  IF len>=0 THEN s[len]:=0C END;  
  lc:=1;
  p:=0;
  ci:=MIN(COLS);
  LOOP
    WHILE (VAL(INTEGER,p)<len) & (s[p]<" ") DO INC(p) END;                    (* remove cr lf *)
    IF (VAL(INTEGER,p)<len) &  (s[p]<>"#") THEN
      card(s, p, colours[ci].r, len, ok);
      IF ok THEN card(s, p, colours[ci].g, len, ok) END;
      IF ok THEN card(s, p, colours[ci].b, len, ok) END;
      IF NOT ok THEN
        IntToStr(lc, 1, h); Append(h, ": colour file syntax error"+LF); 
        Werr(h);
      EXIT END;

      IF ci>=MAX(COLS) THEN EXIT END;

      INC(ci);

    END;
    WHILE (VAL(INTEGER,p)<len) & (s[p]>=" ") DO INC(p) END;                   (* got to eol *)

    INC(lc); 
    IF VAL(INTEGER,p)>=len THEN EXIT END;    
  END; 
END  readcolours;


PROCEDURE Parms;
VAR err:BOOLEAN;
    h:ARRAY[0..1023] OF CHAR;
    posr:POSITION;
    label:CARDINAL;
BEGIN

(* default colours *)
  FILL(ADR(colours), 0C, SIZE(colours));
  colours[cEARTH].r:=200; colours[cEARTH].g:=120; colours[cEARTH].b:=0;
  colours[cHEAVENup].r:=20; colours[cHEAVENup].g:=40; colours[cHEAVENup].b:=180;
  colours[cHEAVENdown].r:=60; colours[cHEAVENdown].g:=110; colours[cHEAVENdown].b:=180;
  colours[cFRESIN].r:=290; colours[cFRESIN].g:=80; colours[cFRESIN].b:=0;
  colours[cFRESOUT].r:=120; colours[cFRESOUT].g:=180; colours[cFRESOUT].b:=20;
  colours[cOPTALT].r:=90; colours[cOPTALT].g:=150; colours[cOPTALT].b:=150;
  colours[cSCALERS].r:=600; colours[cSCALERS].g:=600;colours[cSCALERS].b:=500;
  colours[cMLINES].r:=25; colours[cMLINES].g:=25; colours[cMLINES].b:=25;
  colours[cFRESLINE].r:=130; colours[cFRESLINE].g:=130; colours[cFRESLINE].b:=0;
  colours[cMEADOW].r:=80; colours[cMEADOW].g:=250; colours[cMEADOW].b:=0;
  colours[cTEXT1].r:=800; colours[cTEXT1].g:=700; colours[cTEXT1].b:=100;
  colours[cTEXT2].r:=800; colours[cTEXT2].g:=700; colours[cTEXT2].b:=500;
  colours[cTEXT3].r:=800; colours[cTEXT3].g:=700; colours[cTEXT3].b:=500;
  colours[cTEXTSCALE].r:=700; colours[cTEXTSCALE].g:=700; colours[cTEXTSCALE].b:=700;
  colours[cBRANCH].r:=0; colours[cBRANCH].g:=500; colours[cBRANCH].b:=0;
  colours[cTRUNC].r:=300; colours[cTRUNC].g:=200; colours[cTRUNC].b:=0;
  colours[cWATER].r:=0; colours[cWATER].g:=200; colours[cWATER].b:=300;
  colours[cWOOD].r:=150; colours[cWOOD].g:=120; colours[cWOOD].b:=50;
  colours[cURBAN].r:=400; colours[cURBAN].g:=50; colours[cURBAN].b:=20;

  err:=FALSE;
  label:=0;
  LOOP
    NextArg(h);
    IF h[0]=0C THEN EXIT END;
    IF (h[0]="-") & (h[1]<>0C) & (h[2]=0C) THEN

     IF h[1]="i" THEN
        NextArg(imagefn);
        IF (imagefn[0]=0C) OR (imagefn[0]="-") THEN Error("-i <imagefilename>") END;
 
     ELSIF h[1]="c" THEN
        NextArg(csvfn);
        IF (csvfn[0]=0C) OR (csvfn[0]="-") THEN Error("-c <csv-filename>") END;

     ELSIF h[1]="C" THEN
        NextArg(h);
        IF (h[0]=0C) OR (h[0]="-") THEN Error("-C <colours-filename>") END;
        readcolours(h);

     ELSIF h[1]="x" THEN
        NextArg(h);
        IF NOT StrToInt(h, xsize) OR (xsize<6*20) OR (xsize>=8192) THEN Error("-x <size>") END;
 
     ELSIF h[1]="y" THEN
        NextArg(h);
        IF NOT StrToInt(h, ysize) OR (ysize<10*5) OR (ysize>=8192) THEN Error("-y <size>") END;

     ELSIF h[1]="p" THEN
        NextArg(srtmdir);
        IF (srtmdir[0]=0C) OR (srtmdir[0]="-") THEN Error("-p <srtm-folder-path>") END;

     ELSIF h[1]="a" THEN
        NextArg(h);
        loctopos(posr, h); posa.lat:=posr.lat; posa.long:=posr.long;
        IF NOT posvalid(posa) THEN 
          IF NOT StrToFixl(posa.lat, h) OR (ABS(posa.lat)>=90) THEN Error("-a <lat> <long> or <locator>") END;
          NextArg(h);
          IF NOT StrToFixl(posa.long, h) OR (ABS(posa.long)>180) THEN Error("-a <lat> <long> or <locator>") END;
          posa.lat :=posa.lat*RAD;
          posa.long:=posa.long*RAD;
        END; 
        label:=0;

     ELSIF h[1]="b" THEN
        NextArg(h);
        loctopos(posr, h); posb.lat:=posr.lat; posb.long:=posr.long;
        IF NOT posvalid(posb) THEN
          IF NOT StrToFixl(posb.lat, h) OR (ABS(posb.lat)>=90) THEN Error("-b <lat> <long> or <locator>") END;
          NextArg(h);
          IF NOT StrToFixl(posb.long, h) OR (ABS(posb.long)>180) THEN Error("-b <lat> <long> or <locator>") END;
          posb.lat :=posb.lat*RAD;
          posb.long:=posb.long*RAD;
        END;
        label:=1;
       
     ELSIF h[1]="A" THEN
        NextArg(h);
        IF NOT StrToFixl(anta, h) OR (anta<0) OR (anta>20000) THEN Error("-A <meter> (0..20000)") END;
     ELSIF h[1]="B" THEN
        NextArg(h);
        IF NOT StrToFixl(antb, h) OR (antb<0) OR (antb>20000) THEN Error("-B <meter> (0..20000)") END;

     ELSIF CAP(h[1])="R"  THEN
        opt:=h[1]="R";
        NextArg(h);
        IF NOT StrToFixl(refraction, h) OR (refraction<0.0) OR (refraction>1.0)
        THEN Error("-r <refraction> [0.0..1.0]") END;

     ELSIF h[1]="f" THEN
        NextArg(h);
        IF NOT StrToFixl(mhz, h) OR (mhz<1.0) THEN Error("-f <MHz>") END;

     ELSIF h[1]="F" THEN
        NextArg(h);
        IF h[0]="1" THEN fonttyp:=6 ELSIF h[0]="2" THEN fonttyp:=8 ELSIF h[0]="3" THEN fonttyp:=10
        ELSE Error("-F <font> (1:6x10 2:10x20)") END;

     ELSIF h[1]="g" THEN
        NextArg(h);
        IF NOT StrToFixl(igamma, h) OR (igamma<0.01) OR (igamma>10.0) THEN Error("-g <gamma> [0.1..10]") END;

     ELSIF h[1]="k" THEN
        NextArg(h);
        IF NOT StrToFixl(kernel, h) OR (kernel<0.0) OR (kernel>0.999) THEN Error("-k <factor> [0.0..0.99]") END;

     ELSIF h[1]="M" THEN
        NextArg(h);
        IF NOT StrToFixl(maxdist, h) OR (maxdist<=MINDIST) OR (maxdist>1000.0) THEN Error("-N <km>") END;

     ELSIF h[1]="w" THEN
        NextArg(h);
        IF NOT StrToFixl(treesize, h) OR (treesize<0.0) OR (treesize>100.0) THEN Error("-w <m> [0..100]") END;

     ELSIF h[1]="L" THEN
        NextArg(h);
        IF label=0 THEN Assign(labela, h) ELSIF label=1 THEN Assign(labelb, h) END; 
        INC(label);

     ELSIF h[1]="h" THEN
        WrStrLn("");
        WrStrLn("Geoprofile from Position A to Position B");
        WrStrLn(" -A <m>                            Antenna A over ground [m] (10)");
        WrStrLn(" -a <lat> <long> | [locator]       Position A lat long (degrees) or qth locator");
        WrStrLn(" -B <m>                            Antenna B over ground [m] (10)");
        WrStrLn(" -b <lat> <long> | [locator]       Position B lat long (degrees) or qth locator");
        WrStrLn(" -C <filename>                     Colours File Name (red green blue 0..1023)");
        WrStrLn(" -c <filename>                     csv File Name");
        WrStrLn(" -F <font>                         Font Size (1) 1: 6x10, 2: 8x14, 3: 10x20");
        WrStrLn(" -f <MHz>                          Frequency for Fresnelzone (145)");
        WrStrLn(" -g <gamma>                        Image Gamma 0.1..10 (2.2)");
        WrStrLn(" -h                                this");
        WrStrLn(" -i <filename>                     Image File Name");
        WrStrLn(" -k <factor>                       Kernel diameter of fresnel zone to full diameter (0.6)");
        WrStrLn(" -L <text>                         Label, apply after -a and -b");
        WrStrLn(" -p <pathname>                     folder with /srtm1 /srtm3 /srtm30");
        WrStrLn(" -R <refraction>                   same but earth curvature added to Ground");
        WrStrLn(" -r <refraction>                   0.0(vacuum)..1.0(earth is a disk) (0.25)");
        WrStrLn(" -w <m>                            Wood higth (0) (from 1000 to 2000NN reduced till 0");
        WrStrLn(" -x <size>                         Image size (600)");
        WrStrLn(" -y <size>                         Image size (400)");
        WrStrLn(" use -a and -p Parameter for Altitude as String output");
        WrStrLn("");
        HALT

      ELSE err:=TRUE END;
    ELSE err:=TRUE END;
    IF err THEN EXIT END;
  END;
  IF err THEN Werr(">"); Werr(h); Werr("< use -h"+LF); HALT END;
END Parms;

PROCEDURE wrcsv;
VAR i:INTEGER;
    f:File;
    s,sl:ARRAY[0..1023] OF CHAR;
BEGIN
  f:=OpenWrite(csvfn);
  IF f<0 THEN s:="Cannot write "; Append(s, csvfn); Error(s) END;

  FOR i:=0 TO linksize-1 DO
    WITH path^[i] DO
      IntToStr(framexl+i,1,sl);
      Append(sl, ",");
      FixToStr(pos.lat/RAD, 6, s); Append(sl, s);
      Append(sl, ",");
      FixToStr(pos.long/RAD, 6, s); Append(sl, s);
      Append(sl, ",");
      IntToStr(VAL(INTEGER,optalt+refrm-alt+0.5),1,s); Append(sl, s);
      Append(sl, ",");
      IntToStr(VAL(INTEGER,fresm+0.5),1,s); Append(sl, s);
      Append(sl, LF);
      WrBin(f, sl, Length(sl));       
    END;
  END;

  Close(f);
END wrcsv;


PROCEDURE azimuthl(posa, posb:POSITIONL):REAL;
VAR posar, posbr:POSITION;
BEGIN
  posar.lat:=posa.lat; posar.long:=posa.long;
  posbr.lat:=posb.lat; posbr.long:=posb.long;
  RETURN azimuth(posar, posbr)
END azimuthl;

PROCEDURE fresnel(a, b, lambda:LONGREAL):LONGREAL;
BEGIN IF (lambda<=0.0) OR (a<=0.0) OR (b<=0.0) THEN RETURN 0.0 ELSE RETURN sqrt(lambda*a*b/(a+b)) END; END fresnel;


PROCEDURE calcpath;
VAR res:REAL;
    posr:POSITION;
    posh:POSITIONL;
    x0, y0, z0, x1, y1, z1, dx, dy, dz, k, ks, refrac,
    stepm, a, lambda:LONGREAL;
    i, j, maxsteps, substeps:INTEGER;
    erra, errb:ARRAY[0..99] OF CHAR;
    att:CARD8;
BEGIN
  alta:=getsrtmlong(posa.lat, posa.long, 1, FALSE, res, att, NIL);

  IF alta>=NOALT THEN errorimg("No Altitude Data at Antanna A") END;
  altb:=getsrtmlong(posb.lat, posb.long, 1, FALSE, res, att, NIL);
  IF altb>=NOALT THEN errorimg("No Altitude Data at Antanna B") END;

  alta:=alta+anta;
  altb:=altb+antb;
(*
  IF alta<altb THEN
    a:=alta; alta:=altb; altb:=a;
    posh:=posa; posa:=posb; posb:=posh;
  END; 
*)

  wgs84s(posa.lat, posa.long, alta*0.001, x0,y0,z0);
  wgs84s(posb.lat, posb.long, altb*0.001, x1,y1,z1);
  ks:=1.0/FLOAT(linksize);
  dx:=(x1-x0);
  dy:=(y1-y0);
  dz:=(z1-z0);
  dist:=sqrt(dx*dx + dy*dy + dz*dz);
  stepm:=dist*ks*1000.0;
  maxsteps:=1+TRUNC(stepm*(1.0/30.0));
  erra:="Distance (";
  FixToStr(dist, 3, errb); Append(erra, errb); Append(erra, "km) ");
  
  IF dist<MINDIST THEN
    Append(erra, "too less (>");
    FixToStr(MINDIST, 3, errb); Append(erra, errb); Append(erra, "km)");
    errorimg(erra);
  END;
  IF dist>maxdist THEN
    Append(erra, "too long (<");
    IntToStr(VAL(INTEGER, maxdist), 1, errb); Append(erra, errb); Append(erra, "km)");
    errorimg(erra);
  END;

  refrac:=-(EARTH - sqrt(EARTH*EARTH + dist*dist*0.25))*(4*1000);
  lambda:=300/mhz;
  FOR i:=0 TO linksize-1 DO
    path^[i].alt:=-20000.0;
    substeps:=maxsteps;
    IF i>0 THEN substeps:=1+TRUNC(stepm*4.0/(1.0+ABS(path^[i-1].alt-path^[i-1].optalt))) END;
    IF substeps>maxsteps THEN substeps:=maxsteps END;
    FOR j:=0 TO substeps-1 DO 
      k:=(FLOAT(i) + FLOAT(j)/FLOAT(substeps))*ks;
      WITH path^[i] DO
        wgs84r(x0+dx*k, y0+dy*k, z0+dz*k, pos.lat, pos.long, optalt);
        a:=getsrtmlong(pos.lat, pos.long, TRUNC(stepm), FALSE, resol, attr, NIL);
        IF a>alt THEN alt:=a END;
      END;
    END;
    WITH path^[i] DO
      k:=(FLOAT(i)+1.0)/FLOAT(linksize);
      optalt:=optalt*1000;
      zero:=(k - k*k)*refrac;
      refrm:=zero*refraction;
      fresm:=fresnel(dist*1000*k, dist*1000*(1.0-k), lambda);

      wood:=0.0;
      IF alt>=1.0 THEN
        IF alt>=1000.0 THEN wood:=treesize*(2000.0-alt)/1000.0;
          IF wood<0.0 THEN wood:=0.0 END;
        ELSE wood:=treesize END;
      END;
    END;
  END;
END calcpath;


PROCEDURE fresnelfree(VAR airshadow, woodshadow:LONGREAL); (* guess freeness of fresnel zone *)
VAR i, wd:INTEGER;
    fz, fs, gnd, ws, as:LONGREAL;
BEGIN
  wd:=VAL(INTEGER, FLOAT(linksize)/dist);                              (* steps per km *)
  IF realwood THEN wd:=VAL(INTEGER, FLOAT(linksize)/dist*0.03) END;    (* steps per 30m *)
  airshadow:=0.0;
  woodshadow:=0.0;
  FOR i:=1 TO linksize-1 DO
    WITH path^[i] DO
      fz:=optalt+refrm-fresm*kernel;
      gnd:=alt;
      fs:=fresm*2*kernel;
      IF fs<1.0 THEN fs:=1.0 END;
      as:=(gnd-fz)/fs;
      IF (path^[i].attr=ATTRWOOD) & (i>=wd) & (i<linksize-wd) THEN   (* wood only 1km/30m from endpoints *)
        ws:=(gnd+wood-fz)/fs;
        IF ws>woodshadow THEN woodshadow:=ws END;
      END;
      IF as>airshadow THEN airshadow:=as END;
    END;
  END;
END fresnelfree;

PROCEDURE drawcolon(x:INTEGER; y0, y1:LONGREAL; cr,cg,cb:INTEGER);

  PROCEDURE pix(x,y:INTEGER; lum:REAL);
  BEGIN
    IF lum>0.0 THEN 
      WITH image^[x][y] DO
        INC(r, TRUNC(FLOAT(cr)*lum));
        INC(g, TRUNC(FLOAT(cg)*lum));
        INC(b, TRUNC(FLOAT(cb)*lum));
      END;
    END;
  END pix;


VAR y:INTEGER;
    fy, fx, fx1, fy1, k, yi ,d:LONGREAL;
BEGIN
  k:=y1-y0;
  d:=0.25/sqrt(1.0+sqr(k));
  fx:=0.0;
  REPEAT
    yi:=y0+k*fx;
    y:=TRUNC(yi);
    fy:=yi-FLOAT(y);
    fy1:=1.0-fy;
    fx1:=1.0-fx;
    pix(x,y+1,  fy *fx1);
    pix(x,y,    fy1*fx1);
    pix(x+1,y+1,fy *fx );
    pix(x+1,y,  fy1*fx );
    fx:=fx+d;
  UNTIL fx>1.0;
END drawcolon;


PROCEDURE rightbound(s-:ARRAY OF CHAR):INTEGER;
BEGIN RETURN xsize-VAL(INTEGER, strsize(fonttyp,s))-2 END rightbound;


PROCEDURE Scale125(pix:LONGREAL):INTEGER;
VAR y,d:INTEGER;
BEGIN
  d:=1;
  y:=1;
  LOOP
    IF pix<FLOAT(y*d) THEN RETURN y*d END;
    IF y=1 THEN y:=2 ELSIF y=2 THEN y:=5 ELSE y:=1; d:=d*10 END;
  END;
END Scale125;


PROCEDURE drawimage;

  PROCEDURE drawtree(x:INTEGER; dalt, treesize:LONGREAL);
  VAR i, y, lasttree, branchdense, branches, bty:INTEGER;
      treelum, wty, wtt:LONGREAL;
  BEGIN
    IF treesize<2.0 THEN RETURN END;

    branchdense:=2+VAL(INTEGER, treesize*0.05);
    branches:=VAL(INTEGER, treesize*0.5+sqrt(treesize));
    treelum:=1.5/(log(treesize));

    FOR y:=0 TO VAL(INTEGER,treesize)-1 DO
      bty:=VAL(INTEGER,dalt+treesize)-y;
      IF bty>=ysize THEN bty:=ysize-1 END;
      WITH image^[x+framexl][bty] DO                        (* trunk of tree *)
        INC(r, VAL(INTEGER, FLOAT(colours[cTRUNC].r)*treelum));
        INC(g, VAL(INTEGER, FLOAT(colours[cTRUNC].g)*treelum));
        INC(b, VAL(INTEGER, FLOAT(colours[cTRUNC].b)*treelum));
      END;

      wty:=1.0+FLOAT(y)*(0.2+2.0/treesize);
      i:=1;
      LOOP

        IF (y-i<0) OR (y-i>=branches) THEN EXIT END;

        IF (y-i) MOD branchdense=0 THEN                     (* branch *)
          wtt:=wty-FLOAT(i);                                (* branchs last pixel brightness *)
          IF wtt<0.0 THEN EXIT END;

          IF wtt>1.0 THEN wtt:=1.0 END;
          IF x>i THEN
            WITH image^[x-i+framexl][bty] DO
              INC(r, VAL(INTEGER, FLOAT(colours[cBRANCH].r)*treelum*wtt));
              INC(g, VAL(INTEGER, FLOAT(colours[cBRANCH].g)*treelum*wtt));
              INC(b, VAL(INTEGER, FLOAT(colours[cBRANCH].b)*treelum*wtt));
            END;
          END;
          IF x+i<xsize-framexl THEN
            WITH image^[x+i+framexl][bty] DO
              INC(r, VAL(INTEGER, FLOAT(colours[cBRANCH].r)*treelum*wtt));
              INC(g, VAL(INTEGER, FLOAT(colours[cBRANCH].g)*treelum*wtt));
              INC(b, VAL(INTEGER, FLOAT(colours[cBRANCH].b)*treelum*wtt));
            END;  
          END;
        END;
        INC(i); 
      END;
    END;
  END drawtree;

  
VAR i, x, y, st, n, m, xc, cr, cg, cb, t1, t2, dk, maxi, lasttree:INTEGER;
    min, max, scale, ao, a0, a1, ah, mpp, w, maxw, airshadow, woodshadow:LONGREAL;
    s,ss,sss:ARRAY[0..99] OF CHAR;
    cc:COLS;

  PROCEDURE sc(y:LONGREAL):LONGREAL;
  BEGIN RETURN (y-min)*scale+FLOAT(frameyd) END sc;

  PROCEDURE placetrees(scale:REAL);

    PROCEDURE ai(x:INTEGER):LONGREAL;
    BEGIN RETURN sc(path^[x].alt) END ai;

  VAR i:INTEGER;
      d, lasta, h:LONGREAL;
  BEGIN
    d:=0.0;   
    lasta:=0.0;
    FOR i:=1 TO linksize-2 DO
      IF (path^[i].attr=ATTRWOOD) & (path^[i].wood*scale>4.0) THEN   (* mark wood with priority at begin and end *)
        h:=(2.0*ai(i)-ai(i-1)-ai(i+1))*2.0;                          (* peak *) 
        IF h<0.0 THEN h:=0.0 END;
        IF NOT path^[i-1].treepri
        & (d<=LFLOAT(10*ORD(path^[i+1].attr<>ATTRWOOD)) + h) THEN
          path^[i].treepri:=TRUE;
          h:=path^[i].wood*scale*0.15 (*-ABS(lasta-ai(i))*);
          IF h<0.0 THEN h:=0.0 END;
          d:=d+1.0+h;
          lasta:=ai(i);
        ELSE d:=d-1.0 END;
      ELSE d:=0.0; lasta:=0.0 END;
    END;

  END placetrees;


BEGIN

  fresnelfree(airshadow, woodshadow);
  background;


(*** bottom line ***)
  mpp:=dist*1000.0/VAL(LONGREAL, xsize-1-FRAMEXR-framexl);
  w:=0.0;
  FOR i:=framexl TO xsize-1-FRAMEXR DO
    WITH image^[i][frameyd] DO
      w:=w+mpp;
      IF w>path^[i-framexl].resol THEN
        w:=0.0; 
        r:=600;
        g:=500;
        b:=500;
      ELSE
        r:=300;
        g:=0;
        b:=0;
      END;
    END;
  END;


(*** get scaling ***)
  min:=path^[0].alt;
  max:=min;
  maxw:=min;
  FOR i:=1 TO linksize-1 DO
    WITH path^[i] DO
      optalt:=optalt+refrm;
      IF opt THEN
        alt:=alt+zero;
        optalt:=optalt+zero;
      END;
      IF alt>optalt+fresm THEN wood:=0.0 END;
      IF alt<NOALT THEN
        IF alt<min THEN min:=alt END;
        IF alt>max THEN max:=alt END;
        IF alt+wood>maxw THEN maxw:=alt+wood END;
      END;
      IF optalt+fresm>max THEN max:=optalt+fresm END;
      IF optalt-fresm<min THEN min:=optalt-fresm END;
    END;
  END;
  ah:=max-min;
  w:=maxw-min;                                        (* scale trees to top of image *)
  IF ah<MINM THEN ah:=VAL(LONGREAL,MINM) END;
  w:=FLOAT(ysize-1-frameyd)/w;
  scale:=FLOAT(ysize-(frameyu+frameyd))/ah;
  IF w<scale THEN scale:=w END;                       (* wood to imagesize or rest to top text *)
  min:=min-4.0/scale;                                 (* space for waterblue *)

(*** draw graphs ***)

  placetrees(scale);
  lasttree:=0;
  FOR i:=0 TO linksize-2 DO
    a0:=path^[i].alt;
    a1:=path^[i+1].alt;
    IF (a0<NOALT) & (a1<NOALT) THEN
      ah:=sc(a0);
      IF path^[i].attr=ATTRWATER THEN cc:=cWATER;
      ELSIF path^[i].attr=ATTRWOOD THEN cc:=cWOOD;
      ELSIF path^[i].attr=ATTRURBAN THEN cc:=cURBAN ELSE cc:=cMEADOW END; 
      drawcolon(i+framexl, ah, sc(a1), colours[cc].r, colours[cc].g, colours[cc].b);
      FOR y:=frameyd+1 TO VAL(INTEGER,ah) DO                  (* fill earth *)
        IF (y+4>VAL(INTEGER,ah)) & (path^[i].attr=ATTRWATER) THEN cc:=cWATER ELSE cc:=cEARTH END;
        WITH image^[i+framexl][y] DO
          INC(r, colours[cc].r);
          INC(g, colours[cc].g);
          INC(b, colours[cc].b);  
        END;
      END;


      IF path^[i].treepri THEN
        drawtree(i, sc(path^[i].alt), path^[i].wood*scale);
        treedrawn:=TRUE;
      END;

(*      
      IF (mhz>=30.0) & (a0<path^[i].optalt+path^[i].fresm*kernel)
      & (a0+path^[i].wood>path^[i].optalt-path^[i].fresm*kernel)
      & (path^[i].wood>path^[i].fresm*(0.5*kernel)) THEN       (* tree is in fresnel zone and 1/4 as high *)
      
        w:=path^[i].wood*scale;                                     (* tree size in pixel *)  
        IF (w>4.0) & (i>=lasttree) THEN                             (* space between trees *)
          placetree(lasttree, w);
          IF lasttree<VAL(INTEGER, HIGH(path^)) THEN
            drawtree(lasttree, sc(path^[lasttree].alt), path^[lasttree].wood*scale);
            treedrawn:=TRUE;
          END;
        END;

      END;
*)
    END;
    ao:=path^[i].optalt;
    drawcolon(i+framexl, sc(ao), sc(path^[i+1].optalt), colours[cOPTALT].r, colours[cOPTALT].g, colours[cOPTALT].b);

(*** fresnel as definined ***)
    drawcolon(i+framexl, sc(ao+path^[i].fresm), sc(path^[i+1].optalt+path^[i+1].fresm), colours[cFRESLINE].r, colours[cFRESLINE].g, colours[cFRESLINE].b);
    drawcolon(i+framexl, sc(ao-path^[i].fresm), sc(path^[i+1].optalt-path^[i+1].fresm), colours[cFRESLINE].r, colours[cFRESLINE].g, colours[cFRESLINE].b);
(*
(*** fresnel as effects ***)
    drawcolon(i+framexl, sc(ao+path^[i].fresm*kernel), sc(path^[i+1].optalt+path^[i+1].fresm*kernel), 140, 100, 0);
    drawcolon(i+framexl, sc(ao-path^[i].fresm*kernel), sc(path^[i+1].optalt-path^[i+1].fresm*kernel), 140, 100, 0);
*)
    FOR y:=VAL(INTEGER, sc(ao-path^[i].fresm))+1 TO VAL(INTEGER,sc(ao+path^[i].fresm)) DO
      WITH image^[i+framexl][y] DO                                (* fill fresnel area *)
        IF (y<VAL(INTEGER, sc(ao-path^[i].fresm*kernel))+1)
        OR (y>VAL(INTEGER, sc(ao+path^[i].fresm*kernel))+1) THEN
          INC(r, colours[cFRESOUT].r);
          INC(g, colours[cFRESOUT].g);
          INC(b, colours[cFRESOUT].b);
        ELSE
          INC(r, colours[cFRESIN].r);
          INC(g, colours[cFRESIN].g);
          INC(b, colours[cFRESIN].b);
        END;
      END;
    END;
  END;

(*** vertical scale ***)
  maxi:=linksize-1;
  IF sc(path^[0].optalt)>sc(path^[linksize-1].optalt) THEN maxi:=0 END;
  a0:=path^[maxi].optalt-min;                                       (* meters on scale *)
  a1:=sc(path^[maxi].optalt)-FLOAT(frameyd);                        (* pixel on scale *)
  IF a1=0.0 THEN mpp:=1.0 ELSE mpp:=a0/a1 END;                                                       (* meters per pixel *)
  st:=Scale125(mpp*VAL(LONGREAL,fonty*16+ysize DIV 4)*0.05);  
  m:=(VAL(INTEGER, min)+st-1) DIV st*st;                            (* lowest meters *)

  FOR y:=frameyd TO VAL(INTEGER,sc(path^[maxi].optalt)) DO          (* vertical scale line left *)
    WITH image^[framexl][y] DO
      INC(r, colours[cSCALERS].r);
      INC(g, colours[cSCALERS].g);
      INC(b, colours[cSCALERS].b);
    END;
  END;
  FOR y:=frameyd TO VAL(INTEGER,sc(path^[linksize-1].optalt)) DO    (* vertical scale line right *)
    WITH image^[linksize+(framexl-1)][y] DO
      INC(r, colours[cSCALERS].r);
      INC(g, colours[cSCALERS].g);
      INC(b, colours[cSCALERS].b);
    END;
  END;

  IF opt THEN i:=-3 ELSE i:=0 END;
  LOOP
    y:=frameyd+VAL(INTEGER, (VAL(LONGREAL, st*i+m)-min)/mpp+0.5);
    IF y>=VAL(INTEGER, sc(path^[maxi].optalt)) THEN EXIT END; 

    FOR x:=1 TO linksize-2 DO
      IF opt THEN
        a0:=path^[x].zero*scale+FLOAT(y);
        a1:=path^[x+1].zero*scale+FLOAT(y);
      ELSE a0:=FLOAT(y); a1:=a0 END; 
      IF (a0>FLOAT(frameyd+1)) & (a0 < sc(path^[x].alt)) & (path^[x].alt<NOALT) THEN
        drawcolon(x+framexl, a0, a1, colours[cMLINES].r, colours[cMLINES].g, colours[cMLINES].b);
      END;
    END;

    IF y>frameyd DIV 2+2 THEN
      FOR x:=framexl-2 TO framexl+2 DO                              (* meter lines *)
        WITH image^[x][y] DO
          INC(r, colours[cSCALERS].r);
          INC(g, colours[cSCALERS].g);
          INC(b, colours[cSCALERS].b);
        END;
      END;
      IntToStr(st*i+m, 4, s); Append(s, "m");
      writestr(image, LEFTSPACE, y-fonty DIV 2, fonttyp, 0,0, colours[cTEXTSCALE].r,colours[cTEXTSCALE].g,colours[cTEXTSCALE].b, s);
    END;

    INC(i);
  END;

(*** hor scale ***)
  a0:=dist/FLOAT(linksize)*7000*FLOAT(fontx);
  IF (dist>=10.0) & (a0<1000) THEN a0:=1000.0 END;
  st:=Scale125(a0);
  xc:=VAL(INTEGER, dist*1000) DIV st;
  FOR i:=0 TO xc DO
    n:=framexl+VAL(INTEGER, FLOAT(linksize*i*st)/(dist*1000)+0.5);
    FOR y:=frameyd-2 TO frameyd+1 DO                              (* km lines *)
      WITH image^[n][y] DO
        INC(r, colours[cSCALERS].r);
        INC(g, colours[cSCALERS].g);
        INC(b, colours[cSCALERS].b);
      END;
    END;

    IF st>=1000 THEN dk:=1000 ELSE dk:=1 END;
    IntToStr(st*i DIV dk, 1, s);
    IF i=xc THEN
      IF dk=1 THEN Append(s, "m") ELSE Append(s, "km") END;
    END;
    dk:=n-VAL(INTEGER, strsize(fonttyp,s) DIV 2);
    IF dk+VAL(INTEGER, strsize(fonttyp,s))>=xsize THEN dk:=xsize-VAL(INTEGER, strsize(fonttyp,s))-2 END;
    writestr(image, dk, frameyd-fonty-3, fonttyp, 0,0, colours[cTEXTSCALE].r,colours[cTEXTSCALE].g,colours[cTEXTSCALE].b, s);
  END; 


  IF labela<>"" THEN Assign(s, labela); Append(s, " ["); ELSE Assign(s, "[") END;

  IntToStr(VAL(INTEGER, alta-anta), 0, ss); Append(s, ss); Append(s, "NN+");
  IntToStr(VAL(INTEGER, anta), 0, ss); Append(s, ss); Append(s, "m");
  IF TRUE THEN
    Append(s, " "); FixToStr(azimuthl(posa, posb), 2, ss); Append(s, ss); Append(s, 177C);
  END;
  Append(s, "]");
  writestr(image, LEFTSPACE, ysize-fonty*12 DIV 10, fonttyp, 0,0, colours[cTEXT1].r,colours[cTEXT1].g,colours[cTEXT1].b, s);

  s:="[";
  FixToStr(posa.lat/RAD, 6, ss); Append(s, ss);
  FixToStr(posa.long/RAD, 6, ss); Append(s, "/"); Append(s, ss);
  Append(s, "]");
  writestr(image, LEFTSPACE, ysize-fonty*24 DIV 10, fonttyp, 0,0, colours[cTEXT2].r,colours[cTEXT2].g,colours[cTEXT2].b, s);

  FixToStr(mhz, 2*ORD(mhz<30.0), s); Append(s, "MHz Refrac=");
  FixToStr(refraction, 3, ss); Append(s, ss);
  writestr(image, LEFTSPACE, ysize-fonty*36 DIV 10, fonttyp, 0,0, colours[cTEXT3].r,colours[cTEXT3].g,colours[cTEXT3].b, s);


  Assign(s, "[");
  IntToStr(VAL(INTEGER, altb-antb), 0, ss); Append(s, ss); Append(s, "NN+");
  IntToStr(VAL(INTEGER, antb), 0, ss); Append(s, ss); Append(s, "m");
  IF TRUE THEN
    Append(s, " "); FixToStr(azimuthl(posb, posa), 2, ss); Append(s, ss); Append(s, 177C);
  END;
  Append(s, "]");
  IF labelb<>"" THEN Append(s, " "); Append(s, labelb) END;
  writestr(image, rightbound(s), ysize-fonty*12 DIV 10, fonttyp, 0,0, colours[cTEXT1].r,colours[cTEXT1].g,colours[cTEXT1].b, s);

  s:="[";
  FixToStr(posb.lat/RAD, 6, ss); Append(s, ss);
  FixToStr(posb.long/RAD, 6, ss); Append(s, "/"); Append(s, ss);
  Append(s, "]");
  writestr(image, rightbound(s), ysize-fonty*24 DIV 10, fonttyp, 0,0, colours[cTEXT2].r,colours[cTEXT2].g,colours[cTEXT2].b, s);

--  s:="fspl=";
  FixToStr(32.2+8.68588963806503655305*log(dist*mhz), 2, sss);
  Append(sss,"dB");
  cr:=colours[cTEXT3].r; cg:=colours[cTEXT3].g; cb:=colours[cTEXT3].b;
  IF airshadow>MINAIRSHADOW THEN                                     (* ground shadow*)
    cr:=800; 
    IF airshadow>MAXAIRSHADOW THEN cg:=0;
    ELSE 
      cg:=VAL(INTEGER, (MAXAIRSHADOW-airshadow)*(MAXAIRSHADOW-airshadow)
      *(800/((MAXAIRSHADOW-MINAIRSHADOW)*(MAXAIRSHADOW-MINAIRSHADOW))));
      IF cg<0 THEN cg:=0 ELSIF cg>800 THEN cg:=800 END;
    END;
    cb:=cg;
  END;
  IF (airshadow<0.2) & (woodshadow>0.2) THEN cr:=800; cg:=800; cb:=0 END;             (* wood shadow *)
  writestr(image, rightbound(sss), ysize-fonty*36 DIV 10, fonttyp, 0,0, cr, cg, cb, sss);

  s:="";
  IF treedrawn THEN
    s:="Trees=";
    FixToStr(treesize, 0, ss); Append(s, ss); Append(s, "m ");
  END;
  FixToStr(dist, 2, ss); Append(s, ss); Append(s, "km");
  Append(sss, s); Append(sss, " ");
  writestr(image, rightbound(sss), ysize-fonty*36 DIV 10, fonttyp, 0,0, colours[cTEXT3].r,colours[cTEXT3].g,colours[cTEXT3].b, s);
--WrFixed(airshadow,2, 9);WrFixed(woodshadow,2, 9);WrStrLn(" air wood");

END drawimage;

PROCEDURE wralt;
VAR att:CARD8;
    a, res:REAL;
    
BEGIN
   a:=getsrtmlong(posa.lat, posa.long, 1, FALSE, res, att, NIL);
   IF a>=NOALT THEN Error("no altitude for this position") END;
   WrInt(VAL(INTEGER, a+0.5), 1);
END wralt;


BEGIN
  imagefn:="";
  csvfn:="";
  image:=NIL;
  xsize:=600;
  ysize:=400;
  posinval(posa);
  posinval(posb);
  anta:=10.0;
  antb:=10.0; 
  refraction:=0.25;
  mhz:=145.0;
  opt:=FALSE;
  igamma:=2.2;
  treesize:=0.0;
  fonttyp:=6;
  labela:="";
  labelb:="";
  maxdist:=DEFMAXDIST;
  kernel:=FRESPERC;
  realwood:=TRUE;
  Parms;
  fontx:=fontsizex(fonttyp);
  fonty:=fontsizey(fonttyp);
  framexl:=fontx*5+2+LEFTSPACE;
  frameyd:=fonty+5;
  frameyu:=fonty*39 DIV 10;
  treedrawn:=FALSE;
  IF (posa.lat=0.0) & (posa.long=0.0) THEN Error("need Position A") END;
  IF (posb.lat=0.0) & (posb.long=0.0) THEN
    IF imagefn[0]=0C THEN wralt ELSE Error("need Position B") END;
  ELSE
    IF (imagefn[0]=0C) & (csvfn[0]=0C) THEN Error("need Image and/or csv Filename") END;
    IF xsize<framexl+FRAMEXR+30 THEN Error("xsize too less") END;
    IF ysize<30 THEN Error("ysize too less") END;

    linksize:=xsize-(framexl+FRAMEXR);
    NEW(path, linksize);
    IF path=NIL THEN Error("out of memory") END;  

    calcpath;
    closesrtmfile;                                          (* free srtm cache *) 


    IF csvfn[0]<>0C THEN wrcsv END;

    IF imagefn[0]<>0C THEN
      NEW(image, xsize, ysize);
      IF image=NIL THEN Error("out of memory") END;

      drawimage;

      wrpng;
    END;
  END;
END profile.

(*
gcc  -o profile Lib.o aprspos.o aprsstr.o filesize.o flush.o libsrtm.o osi.o pngwrite.o profile.o /usr/local/xds/lib/x86/libts.a /usr/local/xds/lib/x86/libxds.a -lm libpng.a libz.a imagetext.o
*)