home *** CD-ROM | disk | FTP | other *** search
/ ftp.barnyard.co.uk / 2015.02.ftp.barnyard.co.uk.tar / ftp.barnyard.co.uk / cpm / walnut-creek-CDROM / CPM / HAMRADIO / MICROMUF.PAS < prev    next >
Pascal/Delphi Source File  |  2000-06-30  |  8KB  |  275 lines
  1. PROGRAM MICROMUF (INPUT,OUTPUT);
  2. {Micromuf - A program to computer the minimum and maximum usable frequencies
  3. of a shortwave path between two specified coordinates. }
  4.  
  5. CONST PI = 3.14159265;
  6.        I = 'INVALID';
  7.  
  8. VAR a,n,yt,xt,yr,xr,r,x,h,y,u,q,la,ls,hp,sf,fe,se,re,cp,k,l,xz,mf,ff :REAL;
  9.     xh,z,fo,yf,ty,tl,yg,zo,yz,yn,yi,xn,sx,xs,wx,lh,lm,ab,rd,d,ex,man :REAL;
  10.     transmitter,receiver : string[20];
  11.     m,mh,ve,ho,t :integer;
  12.     correct :boolean;
  13.     key : char;
  14.  
  15. Function power(man,ex:real):real;
  16.    Begin
  17.       power:=EXP(ex*LN(man));
  18.    End;
  19.  
  20. Procedure interlat; { Intermediate Latitude & Longitude calculations }
  21.    Begin
  22.       q:=cos(u*rd)*cos(xt*rd)*sin(k*lm*rd);
  23.       x:=q+sin(xt*rd)*cos(k*lm*rd);
  24.       xn:=(arctan(x/sqrt(-x*x+1+1E-12)))*d;
  25.       q:=(cos(k*lm*rd)-sin(xt*rd)*sin(xn*rd));
  26.       yi:=(-arctan(x/sqrt(-x*x+1+1E-12))+(PI/2))*d;
  27.       if u < 180.0 then yi:=-yi;
  28.       yn:=yt+yi;
  29.       if yn > 180.0 then yn:=yn-360.0;
  30.       if yn <-180.0 then yn:=yn+360.0;
  31.    End;
  32.  
  33. Procedure mini_f2;
  34.    Begin
  35.       yz:=yn;
  36.       If yn<-160.0 then yz:=yn+360.0;
  37.       yg:=(20.0-yz)/50.0;
  38.       zo:=20.0*yg/(1+yg+sqr(yg))+5.0*sqr(1-yg/7.0);
  39.       z:=xn-zo;
  40.       tl:=t-yn/15.0;
  41.       if tl > 24.0 then tl:=tl-24.0;
  42.       if tl < 0.0 then tl:=tl+24.0;
  43.       mh:=m;
  44.       If z <= 0.0 then
  45.          Begin
  46.             z:=-z;
  47.             mh:=mh+6;
  48.          End;
  49.       xh:=cos(30.0*(mh-6.5)*rd);         { 1 week delay on equinoxes }
  50.       sx:=(abs(xh)+xh)/2.0;              { F-layer local summer variance}
  51.       wx:=(abs(xh)-xh)/2.0;              { F-layer local winter variance}
  52.       If z > 77.5 then z:= 77.5;
  53.       ty:=tl;
  54.       If ty < 5.0 then ty:=tl+24.0;
  55.       yf:=(ty-14.0-sx*2.0+wx*2.0-r/175.0)*(7.0-sx*3.0+wx*4.0-r/(150.0-wx*75.0));
  56.       If abs(yf) > 60.0 then yf:=60.0;
  57.       x:=(1+r/(175.0+sx*175.0));
  58.       fo:=6.5*x*cos(yf*rd)*sqrt(cos((z-sx*5.0+wx*5.0)*rd));
  59.       ex:=-0.5;
  60.       sf:=power((1.0-sqr(cos(a*rd)*6367.0/(6367.0+h))),ex);
  61.       ff:=fo*sf;
  62.    End;
  63.  
  64. Procedure e_layer;
  65.    Begin
  66.       q:=sin(xn*rd)*sin(xs*rd);
  67.       x:=q+cos(xn*rd)*cos(xs*rd)*cos((yn-15.0*(t-12.0))*rd);
  68.       xz:=(-arctan(x/sqrt(-x*x+1+1E-12))+PI/2)*d;
  69.       If xz <= 85.0 then
  70.          Begin
  71.             ex:=(1.0/3.0);
  72.             fe:=3.4*(1.0+0.0016*r)*power(cos(xz*rd),ex);
  73.          End
  74.          Else
  75.             Begin
  76.                ex:=-0.5;
  77.                fe:=3.4*(1.0+0.0016*r)*power((xz-80.0),ex);
  78.             End;
  79.       se:=power(1.0-0.965*sqr(cos(a*rd)),ex);
  80.       ls:=0.028*sqr(fe)*se;
  81.    End;
  82.  
  83. Begin { Main Program }
  84.    rd:=PI/180;
  85.    d:=180/PI;
  86.    correct:=FALSE;
  87.    ClrScr;
  88.    Writeln ('                              *** MICROMUF ***     ');
  89.    Writeln;
  90.    Writeln (' This program calculate the :');
  91.    Writeln (' * M. U. F. (Maximum Usable Frequency)');
  92.    Writeln;
  93.    Writeln (' * L. U. F. (Lowest Usable Frequency)');
  94.    Writeln;
  95.    Writeln (' of any shortwave sky-wave path.');
  96.    Writeln;
  97.    Writeln (' Calculations can be done for any month and sunspot number.');
  98.    Writeln;
  99.    Writeln;
  100.    Writeln ('Name transmitter location');
  101.    Readln (transmitter);
  102.    Writeln;
  103.    Repeat;
  104.       Writeln ('Transmitter longitude in degrees. (W=+, E=-)');
  105.       Readln (yt);
  106.       If (yt >=-180.0) and (yt <= 180.0) Then correct:=TRUE
  107.       Else Writeln(I);
  108.    Until correct = TRUE;
  109.    correct:=FALSE;
  110.    Writeln;
  111.    Repeat;
  112.       Writeln ('Transmitter lattitude in degrees. (N=+, S=-)');
  113.       Readln(xt);
  114.       If (xt >= -90.0) and (xt <= 90.0) Then correct:=TRUE
  115.       Else Writeln(I);
  116.    Until correct = TRUE;
  117.    correct:=FALSE;
  118.    Writeln;
  119.    Writeln ('Name receiver location.');
  120.    Readln(receiver);
  121.    Writeln;
  122.    Repeat;
  123.       Writeln ('Receiver longitude in degrees. (W=+, E=-)');
  124.       Readln(yr);
  125.       If (yr >= -180.0) and (yr <=180.0) Then correct:=TRUE
  126.       Else Writeln(I);
  127.    Until correct = TRUE;
  128.    correct:=FALSE;
  129.    Writeln;
  130.    Repeat;
  131.       Writeln ('Receiver lattitude in degrees. (N=+, S=-)');
  132.       readln(xr);
  133.       If (xr >=-90.0) and (xr<=90.0) Then correct:=TRUE
  134.       Else Writeln(I);
  135.    Until correct = TRUE;
  136.    correct:=FALSE;
  137.    Writeln;
  138.    Repeat;
  139.       Writeln ('Sunspot number.');
  140.       Readln (r);
  141.       If (r >= 1.0) and (r <=180.0) Then correct:=TRUE
  142.       Else Writeln(I);
  143.    Until correct = TRUE;
  144.    correct:=FALSE;
  145.    Writeln;
  146.    Repeat;
  147.       Writeln ('Month.');
  148.       Readln (m);
  149.       If (m >= 1) and (m <= 12) Then correct:=TRUE
  150.       Else Writeln(I);
  151.    Until correct = TRUE;
  152.  
  153. {   Geometry Calculations  }
  154.  
  155.    q:=sin(xt*rd)*sin(xr*rd);
  156.    x:=q+cos(xt*rd)*cos(xr*rd)*cos(yt*rd-yr*rd);
  157.    la:=(-arctan(x/sqrt(-x*x+1+1E-12))+(PI/2))*d;
  158.    l:=111.1*la;
  159.    q:=(sin(xr*rd)-sin(xt*rd)*cos(la*rd));
  160.    x:=q/cos(xt*rd)/sin(la*rd);
  161.    u:=(-arctan(x/sqrt(-x*x+1+1E-12))+(PI/2))*d;
  162.    If yt-yr <= 0 Then u:=360-u;
  163.    h:=275+r/2;
  164.    xs:=23.4*cos(30*(m-6.25)*rd);
  165.    n:=n+1;
  166.    lh:=l/n;
  167.    While lh > 4500.0 Do
  168.       Begin
  169.          n:=n+1;
  170.          lh:=l/n;
  171.       End;  { While }
  172.    lm:=la/n;
  173.    a:=(arctan((cos(0.5*lm*rd)-6367.0/(h+6367.0))/sin(0.5*lm*rd)))*d;
  174.    While a < 1.5 Do
  175.       Begin
  176.          n:=n+1;
  177.          lh:=lh/n;
  178.          While lh > 4500.0 Do
  179.             Begin
  180.                n:=n+1;
  181.                lh:=l/n;
  182.             End; { While }
  183.          lm:=la/n;
  184.          a:=(arctan((cos(0.5*lm*rd)-6367.0/(h+6367.0))/sin(0.5*lm*rd)))*d;
  185.       End;
  186.  
  187. { Plot chart on screen }
  188.  
  189.    ClrScr;
  190.    Writeln ('From: ',transmitter,' to: ',receiver);
  191.    Write   ('Month: ',m);
  192.    Writeln (' SSN: ',r:3:0,' Dist: ',Round(l+0.5),' KM');
  193.    Writeln ('Azim: ',Round(u+0.5),' degrees.  F-hops: ',n:4:2);
  194.    ve:=4;
  195.    ho:=1;
  196.    GotoXY(ho,ve);
  197.    q:=34.0;
  198.    While q >=2.0 Do
  199.       Begin
  200.          Writeln ('I                         I',q:2:0);  { 25 spaces }
  201.          q:=q-2.0;
  202.       End; { While }
  203.    Writeln ('---------------------------');   { 27 dashes }
  204.    Writeln (' 0 2 4 6 8 10  14  18  22 H (UTC)');
  205.    Writeln ('      +: MUF  -: LUF');
  206.    ve:=4;
  207.    ho:=32;
  208.    GotoXY(ho,ve);
  209.    Writeln ('mHz');
  210.    For t:=1 to 24 Do
  211.       Begin
  212.          ab:=0.0;
  213.          k:=0.5;
  214.          interlat;
  215.          mini_f2;
  216.          mf:=ff;
  217.          k:=n-0.5;
  218.          interlat;
  219.          mini_f2;
  220.          If ff < mf then mf := ff;
  221.          ve:=21-Round(mf/2.0+0.5);
  222.          ho:=t+1;
  223.          if ve < 4 then ve:=4;
  224.          GotoXY(ho,ve);
  225.          Writeln ('+');
  226.          While k <= n-0.25 Do
  227.             Begin
  228.                interlat;
  229.                e_layer;
  230.                ab:=ab+ls;
  231.                k:=k+0.5;
  232.             End;
  233.          ve:=20-Round(ab+0.5);
  234.          If ve < 4 Then ve:=4;
  235.          If ve > 20 Then ve:=20;
  236.          GotoXY(ho,ve);
  237.          Writeln('-');
  238.       End;
  239.       Writeln;
  240.       Writeln;
  241.       Writeln;
  242.       Writeln;
  243.       Writeln;
  244.       Write ('          ----- Press a key to continue ----- ');
  245.       Readln (key);
  246. End.
  247.  
  248. {  This program uses 'MINI-F2' devised by R. Fricker (BBC external
  249.   services) for FO-F calculations and L.M. Muggleton's formula for
  250.   FO-E calculations.
  251.  
  252.    For the L.U.F. a minimum useable fieldstrength of 30 DBUV at the
  253.   receiver and 250 KW of transmitter power (aerial gain: 18 DBI) are
  254.   assumed.  The L.U.F. is derived from absorption calculations based
  255.   on the work of Piggot, George, Samuel, and Bradley.  In spite of the
  256.   program's simplicity it gives a good impression of the ionosphere's
  257.   behaviour and can be used for propagation predictions.
  258.  
  259.             Hans Bakhuizen
  260.             Propagation Unit; Frequency Bureau
  261.             Radio Netherlands
  262.             P.O. Box 222
  263.             1200 JG Hilversum Holland
  264.  
  265. (C) Copyright Media Network June 1984
  266.  
  267. Translation from Basic into TURBO Pascal by Jonathan D Ogden on
  268. September 26, 1986.
  269.  
  270.        Jonathan D Ogden
  271.        402 e Daniel
  272.        Champaign, Il 61820 USA
  273.  
  274. }
  275.