<p>There were a few other small errors in your translation (such as indenting too much stuff in two of the <code>if</code> statements), and I think you might have been trying to be too clever with changing <code>SUBROUTINE</code> into defined functions when it was a lot easier to just drop the code in-line since it was only used once.</p>

<p>There were minor tweaks and as far as I can tell this seems to work (you'll want to play with the formatted print to clean it up, and this runs under Python 3; if you need it to run under Python 2, you'll probably need to change the <code>print</code> and <code>num1.write</code> statements):</p>

<pre><code>#!/usr/bin/env python

import math

print('Program to flag Jovian Decametric windows')
month = ['JAN','FEB','MAR','APR','MAY','JUN','JUL','AUG','SEP','OCT','NOV','DEC']
week = 42.46/360
pi = 3.141593
kr = pi / 180
num1 = open('jovrad.txt', 'w')
yy = int(input(("Year for which predictions are required ")))
e = math.trunc(((yy-1) / 100))
print(e)
f = 2 - e + math.trunc(e/4)
print(f)
jd = math.trunc(365.25 * (yy - 1)) + 1721423 + f + .5
print(jd)
d0 = jd - 2435108
print(d0)
incr = 0
dmax = 0
tx = 0
ty = 0
yyly = 0
yylc = 0
if yy / 400 - math.trunc((yy / 400)) == 0:
    incr = 1
    print("in if-1")
yyly = yy / 4 - math.trunc((yy / 4))
yylc = yy / 100 - math.trunc((yy / 100))
if yyly == 0 and yylc != 0:
    print("in if-2")
    incr = 1
ty = 59 + incr
dmax = 365 + incr
tx = ty + .5
num1.write("******************************************************\n")
pout = "  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR " + str(yy) + "\n"
num1.write(pout)
num1.write("******************************************************\n")
num1.write("\n")
num1.write("Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source")
num1.write("\n")
th = 0

while int(th / 24) + 1 &lt;= dmax:
    d = d0 + (th / 24)
    v = (157.0456 + .0011159 * d) % 360
    m = (357.2148 + .9856003 * d) % 360
    n = (94.3455 + .0830853 * d + .33 * math.sin(kr * v)) % 360
    j = (351.4266 + .9025179 * d - .33 * math.sin(kr * v)) % 360
    a = 1.916 * math.sin(kr * m) + .02 * math.sin(kr * 2 * m)
    b = 5.552 * math.sin(kr * n) + .167 * math.sin(kr * 2 * n)
    k = j + a - b
    r = 1.00014 - .01672 * math.cos(kr * m) - .00014 * math.cos(kr * 2 * m)
    re = 5.20867 - .25192 * math.cos(kr * n) - .0061 * math.cos(kr * 2 * n)
    dt = math.sqrt(re * re + r * r - 2 * re * r * math.cos(kr * k))
    sp = r * math.sin(kr * k) / dt
    ps = sp / .017452
    dl = d - dt / 173
    pb = ps - b
    xi = 150.4529 * math.trunc((dl)) + 870.4529 * (dl - math.trunc((dl)))
    L3 = (274.319 + pb + xi + .01016 * 51) % 360
    U1 = 101.5265 + 203.405863 * dl + pb
    U2 = 67.81114 + 101.291632 * dl + pb
    z = (2 * (U1 - U2)) % 360
    U1 = U1 + .472 * math.sin(kr * z)
    U1 = (U1 + 180) % 360
    s = ""
    if L3 &lt; 255 and L3 &gt; 200 and U1 &lt; 250 and U1 &gt; 220:
        s = "Io-A"
    if L3 &lt; 180 and L3 &gt; 105 and U1 &lt; 100 and U1 &gt; 80:
        s = "Io-B"
    if L3 &lt; 350 and L3 &gt; 300 and U1 &lt; 250 and U1 &gt; 230:
        s = "Io-C"
    if s != "":
        dy = math.trunc((th / 24)) + 1
        h = th - (dy - 1) * 24
        if(dy &gt; th):
            m = math.trunc((dy - tx) / 30.6) + 3
            da = dy - ty - math.trunc((m - 3) * 30.6 + .5)
        else:
            m = math.trunc((dy - 1) / 31) + 1
            da = dy - (m - 1) * 31
        mn = month[(m-1)]
#        mn = month[(m-1)*3+1:(m-1)*3-1+3]
        outstring = "%i  %s  %i  %5.3f  %5.3f  %5.3f  %5.3f  %s\n" % (dy, mn, da, h, U1, L3, dt, s)
        num1.write(outstring)
    th = th + .5

num1.close()
print("Program Complete  - results in file JOVRAD.TXT")
</code></pre>

CC BY-SA 3.0
2015-11-25T00:59:35.010Z
Thank you so much! I really appreciate it. This does run 100% in Python 2 by the way. And yeah I'll have to work a little bit on the printing format.
Macuser
CC BY-SA 3.0
2015-11-25T01:07:26.000Z
Also, I removed the first variable in the "outstring" you wrote. It just starts with the month and then the day (variable "da").
Macuser
CC BY-SA 3.0
2015-11-24T23:17:12.950Z
2015-11-24T23:17:12.950Z
Dave
<p>A few days ago I decided to take on the little project of converting a Qbasic program into Python (as a side project to my Radio Jove project), and I've managed to get it to run, but the math is definitely off. I was hoping for a fellow astronomer or astrophysicist to assist in the mathematical side of my Python code. Obviously something is wrong with the math aspect of the code since most of the values in the python program are very skewed. My apologies in advance for the single character variables. I would make them better if I knew everything going on.</p>

<p>Qbasic code:</p>

<pre><code>'Program to flag Jovian Decametric windows
 m$ = "JANFEBMARAPRMAYJUNJULAUGSEPOCTNOVDEC"
     wi = 42.46 / 360
     pi = 3.141593
     kr = pi / 180
     f$ = "###  \ \ ##   ##.#     ###      ###    ##.##     \  \"
 OPEN "jovrad.txt" FOR OUTPUT AS #1
 INPUT "Year for which predictions are required"; yy
 e = INT((yy - 1) / 100)
 f = 2 - e + INT(e / 4)
 jd = INT(365.25 * (yy - 1)) + 1721423 + f + .5
 d0 = jd - 2435108
 incr = 0
 IF yy / 400 - INT(yy / 400) = 0 THEN incr = 1
 yyly = yy / 4 - INT(yy / 4)
 yylc = yy / 100 - INT(yy / 100)
 IF yyly = 0 AND yylc &lt;&gt; 0 THEN incr = 1
 ty = 59 + incr
 dmax = 365 + incr
 tx = ty + .5
 PRINT #1, "******************************************************"
 PRINT #1, "  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR"; yy
 PRINT #1, "******************************************************"
 PRINT #1,
 PRINT #1, "Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source"
 PRINT #1,
 th = 0
 DO
   GOSUB Compute
   s$ = ""
       IF L3 &lt; 255 AND L3 &gt; 200 AND U1 &lt; 250 AND U1 &gt; 220 THEN s$ = "Io-A"
   IF L3 &lt; 180 AND L3 &gt; 105 AND U1 &lt; 100 AND U1 &gt; 80 THEN s$ = "Io-B"
       IF L3 &lt; 350 AND L3 &gt; 300 AND U1 &lt; 250 AND U1 &gt; 230 THEN s$ = "Io-C"
   IF s$ &lt;&gt; "" THEN GOSUB Outdat
   th = th + .5
 LOOP UNTIL INT(th / 24) + 1 &gt; dmax
 PRINT "Program completed - results in file JOVRAD.TXT"
 END

Compute:
 d = d0 + th / 24
 v = (157.0456 + .0011159# * d) MOD 360
 m = (357.2148 + .9856003# * d) MOD 360
 n = (94.3455 + .0830853# * d + .33 * SIN(kr * v)) MOD 360
 j = (351.4266 + .9025179# * d - .33 * SIN(kr * v)) MOD 360
 a = 1.916 * SIN(kr * m) + .02 * SIN(kr * 2 * m)
 b = 5.552 * SIN(kr * n) + .167 * SIN(kr * 2 * n)
 k = j + a - b
 r = 1.00014 - .01672 * COS(kr * m) - .00014 * COS(kr * 2 * m)
 re = 5.20867 - .25192 * COS(kr * n) - .0061 * COS(kr * 2 * n)
 dt = SQR(re * re + r * r - 2 * re * r * COS(kr * k))
 sp = r * SIN(kr * k) / dt
 ps = sp / .017452
 dl = d - dt / 173
 pb = ps - b
 xi = 150.4529 * INT(dl) + 870.4529 * (dl - INT(dl))
 L3 = (274.319 + pb + xi + .01016 * 51) MOD 360
 U1 = 101.5265 + 203.405863# * dl + pb
 U2 = 67.81114 + 101.291632# * dl + pb
 z = (2 * (U1 - U2)) MOD 360
 U1 = U1 + .472 * SIN(kr * z)
 U1 = (U1 + 180) MOD 360
RETURN

Outdat:
 dy = INT(th / 24) + 1
 h = th - (dy - 1) * 24
 IF dy &gt; ty THEN
   m = INT((dy - tx) / 30.6) + 3
   da = dy - ty - INT((m - 3) * 30.6 + .5)
 ELSE
   m = INT((dy - 1) / 31) + 1
   da = dy - (m - 1) * 31
 END IF
 mn$ = MID$(m$, (m - 1) * 3 + 1, 3)
     PRINT #1, USING f$; dy; mn$; da; h; U1; L3; dt; s$
RETURN
</code></pre>

<p>Sample output:</p>

<pre><code>******************************************************
  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR 1994 
******************************************************

Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source

  4  JAN  4   11.0     228      205     5.81     Io-A
  4  JAN  4   11.5     233      224     5.81     Io-A
  4  JAN  4   12.0     237      242     5.81     Io-A
  6  JAN  6    6.0     233      325     5.78     Io-C
  6  JAN  6    6.5     237      343     5.78     Io-C
  7  JAN  7    6.5      81      133     5.76     Io-B
  7  JAN  7    7.0      86      152     5.76     Io-B
  7  JAN  7    7.5      90      170     5.76     Io-B
  9  JAN  9   20.0     241      203     5.73     Io-A
  9  JAN  9   20.5     246      222     5.73     Io-A
 11  JAN 11   12.5     225      232     5.70     Io-A
 11  JAN 11   13.0     229      251     5.70     Io-A
 11  JAN 11   14.5     242      305     5.70     Io-C
 11  JAN 11   15.0     246      323     5.70     Io-C
 12  JAN 12   15.0      90      113     5.68     Io-B
 12  JAN 12   15.5      94      131     5.68     Io-B
 12  JAN 12   16.0      98      150     5.68     Io-B
 14  JAN 14    8.5      82      179     5.66     Io-B
 16  JAN 16   21.0     234      214     5.61     Io-A
 16  JAN 16   21.5     238      232     5.61     Io-A
 16  JAN 16   22.0     242      250     5.61     Io-A
 18  JAN 18   15.5     234      315     5.59     Io-C
 18  JAN 18   16.0     238      333     5.59     Io-C
 19  JAN 19   16.0      82      124     5.57     Io-B
 19  JAN 19   16.5      86      142     5.57     Io-B
 19  JAN 19   17.0      91      160     5.57     Io-B
 19  JAN 19   17.5      95      178     5.57     Io-B
</code></pre>

<p>Python code:</p>

<pre><code>#!/usr/bin/env python
from __future__ import print_function, division
import math
print('Program to flag Jovian Decametric windows')
month = ['JAN','FEB','MAR','APR','MAY','JUN','JUL','AUG','SEP','OCT','NOV','DEC']
week = 42.46/360
pi = 3.141593
kr = pi / 180
form = "###  \ \ ##   ##.#     ###      ###    ##.##     \  \\"
num1 = open('jovrad.txt', 'w')
yy = int(raw_input(("Year for which predictions are required ")))
e = math.trunc(((yy-1) / 100))
f = 2 - e + math.trunc(e/4)
jd = math.trunc(365.25 * (yy - 1)) + 1721423 + f + .5
d0 = jd - 2435108
incr = 0
dmax = 0
tx = 0
ty = 0
yyly = 0
yylc = 0
if yy / 400 - math.trunc((yy / 400)) == 0:
    incr = 1
    yyly = yy / 4 - math.trunc((yy / 4))
    yylc = yy / 100 - math.trunc((yy / 100))
if yyly == 0 and yylc != 0:
    incr = 1
    ty = 59 + incr
    dmax = 365 + incr
    tx = ty + .5
print("******************************************************", file=num1)
print("  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR ",yy, file=num1)
print("******************************************************", file=num1)
print("\n", file=num1)
print("Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source", file=num1)
print("\n", file=num1)
th = 0
def compute(d0, th, dmax):
    global U1, L3, dt, s
    d = d0 + math.trunc(th / 24)
    v = (157.0456 + .0011159 * d) % 360
    m = (357.2148 + .9856003 * d) % 360
    n = (94.3455 + .0830853 * d + .33 * math.sin(kr * v)) % 360
    j = (351.4266 + .9025179 * d - .33 * math.sin(kr * v)) % 360
    a = 1.916 * math.sin(kr * m) + .02 * math.sin(kr * 2 * m)
    b = 5.552 * math.sin(kr * n) + .167 * math.sin(kr * 2 * n)
    k = j + a - b
    r = 1.00014 - .01672 * math.cos(kr * m) - .00014 * math.cos(kr * 2 * m)
    re = 5.20867 - .25192 * math.cos(kr * n) - .0061 * math.cos(kr * 2 * n)
    dt = math.sqrt(re * re + r * r - 2 * re * r * math.cos(kr * k))
    sp = r * math.sin(kr * k) / dt
    ps = sp / .017452
    dl = d - dt / 173
    pb = ps - b
    xi = 150.4529 * math.trunc((dl)) + 870.4529 * (dl - math.trunc((dl)))
    L3 = (274.319 + pb + xi + .01016 * 51) % 360
    U1 = 101.5265 + 203.405863 * dl + pb
    U2 = 67.81114 + 101.291632 * dl + pb
    z = (2 * (U1 - U2)) % 360
    U1 = U1 + .472 * math.sin(kr * z)
    U1 = (U1 + 180) % 360
    s = ""
    while math.trunc(th / 24) + 1 &lt; dmax and math.trunc(th / 24) + 1 == dmax:
        if L3 &lt; 255 and L3 &gt; 200 and U1 &lt; 250 and U1 &gt; 220:
            s = "Io-A"
        if L3 &lt; 180 and L3 &gt; 105 and U1 &lt; 100 and U1 &gt; 80:
            s = "Io-B"
        if L3 &lt; 350 and L3 &gt; 300 and U1 &lt; 250 and U1 &gt; 230:
            s = "Io-C"
        if s != "":
            outdat()
        th = th + .5
    print("Program completed - results in file JOVRAD.TXT", file=num1)
compute(d0, th, dmax)

def outdat(th,tx,ty):
    dy = math.trunc((th / 24)) + 1
    h = th - (dy - 1) * 24
    if(dy &gt; th):
        m = math.trunc((dy - tx) / 30.6) + 3
        da = dy - ty - math.trunc((m - 3) * 30.6 + .5)
    else:
        m = math.trunc((dy - 1) / 31) + 1
        da = dy - (m - 1) * 31
    mn = month[(m-1)*3+1:(m-1)*3-1+3]
    print(dy, mn, da, h, U1, L3, dt, s, file=num1)
outdat(th,tx,ty)
</code></pre>

<p>My current output:</p>

<pre><code>******************************************************
  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR  1998
******************************************************


Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source


-12.9775 1.0 ['AUG'] -0.5 0.0 135.325782651 10.6054462674 5.7071322535 
</code></pre>

CC BY-SA 3.0
2015-11-23T03:09:08.623Z
Nice that you chosed Python. Have you considered using `numpy`? Then you will have all the math functions, like `np.pi`, `np.sin()`, and much more.
skytux
CC BY-SA 3.0
2015-11-23T03:11:22.337Z
Yeah, I really enjoy Python. I do have numpy installed as well as scipy, but I'm not sure if I will need to use it. The Python math library already has a pi function (which I should use in this program) and a sin and cosine function that I use already.
Macuser
CC BY-SA 3.0
2015-11-23T03:13:04.910Z
Ok, but what is your question then? What is wrong with your math? Seems fine to me.
skytux
CC BY-SA 3.0
2015-11-23T03:28:56.500Z
Well considering I got a negative number for the "Day" and that all of the other values are very different then what they should be, something is wrong.
Macuser
CC BY-SA 3.0
2015-11-23T04:13:56.187Z
My bad, I didn't see the output!
skytux
CC BY-SA 3.0
2015-11-23T04:20:17.300Z
@Aabaakawad Yes, the second block of code shows that.
Macuser
CC BY-SA 3.0
2015-11-23T04:33:17.900Z
@Macuser sorry, I'll delete my comment.
Eubie Drew
CC BY-SA 3.0
2015-11-23T11:53:01.597Z
You have explicit casting of results (that are floats) as ints in the basic but not in the Python. I would worry about the differences that result from this. In the Python what types are the variables and how is any implicit casting (if any) to ints done?
Conrad Turner
CC BY-SA 3.0
2015-11-23T23:00:43.357Z
@Contrad Turner is correct.  Take, for instance, the line that assigns something to the variable $e$.  In the QBasic code, you are using the $INT()$ function which truncates the decimal part off of numbers and returns the integer.  To do that in python, you could use the $math.trunc()$ function.  Looking at your code, I think the $INT$ is the only thing you need to take care of.  I believe the $MOD$ works the same as the python modulo, if I recall correctly.
Dave
CC BY-SA 3.0
2015-11-23T23:52:29.930Z
Great! I've updated the python code to have all of the Ints where needed. I ran it again and it seems I'm still getting a negative number for the day value. But I'm sure it fixed some issues that were present before.
Macuser
CC BY-SA 3.0
2015-11-24T02:44:18.593Z
Just updated the code even further. I was accidentally printing out "f" when "dy" was supposed to be the printed as the first thing. Now I just need to fix the formatting so that I can use my "form" variable, and then double check output values.
Macuser
CC BY-SA 3.0
2015-11-24T03:01:36.907Z
The numbers still definitely aren't calculating right. I'm only getting one line of data for any year I enter. Also, the source is still not showing up.
Macuser
<jupiter><radio-astronomy><python><radio>
Need help with the math in Python program to flag Jovian radio emissions
Macuser
# Need help with the math in Python program to flag Jovian radio emissions

(Asked by: username_0 on 2015-11-23)

A few days ago I decided to take on the little project of converting a Qbasic program into Python (as a side project to my Radio Jove project), and I've managed to get it to run, but the math is definitely off. I was hoping for a fellow astronomer or astrophysicist to assist in the mathematical side of my Python code. Obviously something is wrong with the math aspect of the code since most of the values in the python program are very skewed. My apologies in advance for the single character variables. I would make them better if I knew everything going on.
Qbasic code:
'Program to flag Jovian Decametric windows
 m$ = "JANFEBMARAPRMAYJUNJULAUGSEPOCTNOVDEC"
     wi = 42.46 / 360
     pi = 3.141593
     kr = pi / 180
     f$ = "###  \ \ ##   ##.#     ###      ###    ##.##     \  \"
 OPEN "jovrad.txt" FOR OUTPUT AS #1
 INPUT "Year for which predictions are required"; yy
 e = INT((yy - 1) / 100)
 f = 2 - e + INT(e / 4)
 jd = INT(365.25 * (yy - 1)) + 1721423 + f + .5
 d0 = jd - 2435108
 incr = 0
 IF yy / 400 - INT(yy / 400) = 0 THEN incr = 1
 yyly = yy / 4 - INT(yy / 4)
 yylc = yy / 100 - INT(yy / 100)
 IF yyly = 0 AND yylc <> 0 THEN incr = 1
 ty = 59 + incr
 dmax = 365 + incr
 tx = ty + .5
 PRINT #1, "******************************************************"
 PRINT #1, "  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR"; yy
 PRINT #1, "******************************************************"
 PRINT #1,
 PRINT #1, "Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source"
 PRINT #1,
 th = 0
 DO
   GOSUB Compute
   s$ = ""
       IF L3 < 255 AND L3 > 200 AND U1 < 250 AND U1 > 220 THEN s$ = "Io-A"
   IF L3 < 180 AND L3 > 105 AND U1 < 100 AND U1 > 80 THEN s$ = "Io-B"
       IF L3 < 350 AND L3 > 300 AND U1 < 250 AND U1 > 230 THEN s$ = "Io-C"
   IF s$ <> "" THEN GOSUB Outdat
   th = th + .5
 LOOP UNTIL INT(th / 24) + 1 > dmax
 PRINT "Program completed - results in file JOVRAD.TXT"
 END

Compute:
 d = d0 + th / 24
 v = (157.0456 + .0011159# * d) MOD 360
 m = (357.2148 + .9856003# * d) MOD 360
 n = (94.3455 + .0830853# * d + .33 * SIN(kr * v)) MOD 360
 j = (351.4266 + .9025179# * d - .33 * SIN(kr * v)) MOD 360
 a = 1.916 * SIN(kr * m) + .02 * SIN(kr * 2 * m)
 b = 5.552 * SIN(kr * n) + .167 * SIN(kr * 2 * n)
 k = j + a - b
 r = 1.00014 - .01672 * COS(kr * m) - .00014 * COS(kr * 2 * m)
 re = 5.20867 - .25192 * COS(kr * n) - .0061 * COS(kr * 2 * n)
 dt = SQR(re * re + r * r - 2 * re * r * COS(kr * k))
 sp = r * SIN(kr * k) / dt
 ps = sp / .017452
 dl = d - dt / 173
 pb = ps - b
 xi = 150.4529 * INT(dl) + 870.4529 * (dl - INT(dl))
 L3 = (274.319 + pb + xi + .01016 * 51) MOD 360
 U1 = 101.5265 + 203.405863# * dl + pb
 U2 = 67.81114 + 101.291632# * dl + pb
 z = (2 * (U1 - U2)) MOD 360
 U1 = U1 + .472 * SIN(kr * z)
 U1 = (U1 + 180) MOD 360
RETURN

Outdat:
 dy = INT(th / 24) + 1
 h = th - (dy - 1) * 24
 IF dy > ty THEN
   m = INT((dy - tx) / 30.6) + 3
   da = dy - ty - INT((m - 3) * 30.6 + .5)
 ELSE
   m = INT((dy - 1) / 31) + 1
   da = dy - (m - 1) * 31
 END IF
 mn$ = MID$(m$, (m - 1) * 3 + 1, 3)
     PRINT #1, USING f$; dy; mn$; da; h; U1; L3; dt; s$
RETURN

Sample output:
******************************************************
  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR 1994 
******************************************************

Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source

  4  JAN  4   11.0     228      205     5.81     Io-A
  4  JAN  4   11.5     233      224     5.81     Io-A
  4  JAN  4   12.0     237      242     5.81     Io-A
  6  JAN  6    6.0     233      325     5.78     Io-C
  6  JAN  6    6.5     237      343     5.78     Io-C
  7  JAN  7    6.5      81      133     5.76     Io-B
  7  JAN  7    7.0      86      152     5.76     Io-B
  7  JAN  7    7.5      90      170     5.76     Io-B
  9  JAN  9   20.0     241      203     5.73     Io-A
  9  JAN  9   20.5     246      222     5.73     Io-A
 11  JAN 11   12.5     225      232     5.70     Io-A
 11  JAN 11   13.0     229      251     5.70     Io-A
 11  JAN 11   14.5     242      305     5.70     Io-C
 11  JAN 11   15.0     246      323     5.70     Io-C
 12  JAN 12   15.0      90      113     5.68     Io-B
 12  JAN 12   15.5      94      131     5.68     Io-B
 12  JAN 12   16.0      98      150     5.68     Io-B
 14  JAN 14    8.5      82      179     5.66     Io-B
 16  JAN 16   21.0     234      214     5.61     Io-A
 16  JAN 16   21.5     238      232     5.61     Io-A
 16  JAN 16   22.0     242      250     5.61     Io-A
 18  JAN 18   15.5     234      315     5.59     Io-C
 18  JAN 18   16.0     238      333     5.59     Io-C
 19  JAN 19   16.0      82      124     5.57     Io-B
 19  JAN 19   16.5      86      142     5.57     Io-B
 19  JAN 19   17.0      91      160     5.57     Io-B
 19  JAN 19   17.5      95      178     5.57     Io-B

Python code:
#!/usr/bin/env python
from __future__ import print_function, division
import math
print('Program to flag Jovian Decametric windows')
month = ['JAN','FEB','MAR','APR','MAY','JUN','JUL','AUG','SEP','OCT','NOV','DEC']
week = 42.46/360
pi = 3.141593
kr = pi / 180
form = "###  \ \ ##   ##.#     ###      ###    ##.##     \  \\"
num1 = open('jovrad.txt', 'w')
yy = int(raw_input(("Year for which predictions are required ")))
e = math.trunc(((yy-1) / 100))
f = 2 - e + math.trunc(e/4)
jd = math.trunc(365.25 * (yy - 1)) + 1721423 + f + .5
d0 = jd - 2435108
incr = 0
dmax = 0
tx = 0
ty = 0
yyly = 0
yylc = 0
if yy / 400 - math.trunc((yy / 400)) == 0:
    incr = 1
    yyly = yy / 4 - math.trunc((yy / 4))
    yylc = yy / 100 - math.trunc((yy / 100))
if yyly == 0 and yylc != 0:
    incr = 1
    ty = 59 + incr
    dmax = 365 + incr
    tx = ty + .5
print("******************************************************", file=num1)
print("  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR ",yy, file=num1)
print("******************************************************", file=num1)
print("\n", file=num1)
print("Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source", file=num1)
print("\n", file=num1)
th = 0
def compute(d0, th, dmax):
    global U1, L3, dt, s
    d = d0 + math.trunc(th / 24)
    v = (157.0456 + .0011159 * d) % 360
    m = (357.2148 + .9856003 * d) % 360
    n = (94.3455 + .0830853 * d + .33 * math.sin(kr * v)) % 360
    j = (351.4266 + .9025179 * d - .33 * math.sin(kr * v)) % 360
    a = 1.916 * math.sin(kr * m) + .02 * math.sin(kr * 2 * m)
    b = 5.552 * math.sin(kr * n) + .167 * math.sin(kr * 2 * n)
    k = j + a - b
    r = 1.00014 - .01672 * math.cos(kr * m) - .00014 * math.cos(kr * 2 * m)
    re = 5.20867 - .25192 * math.cos(kr * n) - .0061 * math.cos(kr * 2 * n)
    dt = math.sqrt(re * re + r * r - 2 * re * r * math.cos(kr * k))
    sp = r * math.sin(kr * k) / dt
    ps = sp / .017452
    dl = d - dt / 173
    pb = ps - b
    xi = 150.4529 * math.trunc((dl)) + 870.4529 * (dl - math.trunc((dl)))
    L3 = (274.319 + pb + xi + .01016 * 51) % 360
    U1 = 101.5265 + 203.405863 * dl + pb
    U2 = 67.81114 + 101.291632 * dl + pb
    z = (2 * (U1 - U2)) % 360
    U1 = U1 + .472 * math.sin(kr * z)
    U1 = (U1 + 180) % 360
    s = ""
    while math.trunc(th / 24) + 1 < dmax and math.trunc(th / 24) + 1 == dmax:
        if L3 < 255 and L3 > 200 and U1 < 250 and U1 > 220:
            s = "Io-A"
        if L3 < 180 and L3 > 105 and U1 < 100 and U1 > 80:
            s = "Io-B"
        if L3 < 350 and L3 > 300 and U1 < 250 and U1 > 230:
            s = "Io-C"
        if s != "":
            outdat()
        th = th + .5
    print("Program completed - results in file JOVRAD.TXT", file=num1)
compute(d0, th, dmax)

def outdat(th,tx,ty):
    dy = math.trunc((th / 24)) + 1
    h = th - (dy - 1) * 24
    if(dy > th):
        m = math.trunc((dy - tx) / 30.6) + 3
        da = dy - ty - math.trunc((m - 3) * 30.6 + .5)
    else:
        m = math.trunc((dy - 1) / 31) + 1
        da = dy - (m - 1) * 31
    mn = month[(m-1)*3+1:(m-1)*3-1+3]
    print(dy, mn, da, h, U1, L3, dt, s, file=num1)
outdat(th,tx,ty)

My current output:
******************************************************
  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR  1998
******************************************************


Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source


-12.9775 1.0 ['AUG'] -0.5 0.0 135.325782651 10.6054462674 5.7071322535

**Tags:** jupiter, radio-astronomy, python, radio

> **Comments:**
> **username_1**: Nice that you chosed Python. Have you considered using `numpy`? Then you will have all the math functions, like `np.pi`, `np.sin()`, and much more.  
> **username_0**: Yeah, I really enjoy Python. I do have numpy installed as well as scipy, but I'm not sure if I will need to use it. The Python math library already has a pi function (which I should use in this program) and a sin and cosine function that I use already.  
> **username_1**: Ok, but what is your question then? What is wrong with your math? Seems fine to me.  
> **username_0**: Well considering I got a negative number for the "Day" and that all of the other values are very different then what they should be, something is wrong.  
> **username_1**: My bad, I didn't see the output!  
> **username_0**: @Aabaakawad Yes, the second block of code shows that.  
> **username_2**: @username_0 sorry, I'll delete my comment.  
> **username_3**: You have explicit casting of results (that are floats) as ints in the basic but not in the Python. I would worry about the differences that result from this. In the Python what types are the variables and how is any implicit casting (if any) to ints done?  
> **username_4**: @Contrad Turner is correct.  Take, for instance, the line that assigns something to the variable $e$.  In the QBasic code, you are using the $INT()$ function which truncates the decimal part off of numbers and returns the integer.  To do that in python, you could use the $math.trunc()$ function.  Looking at your code, I think the $INT$ is the only thing you need to take care of.  I believe the $MOD$ works the same as the python modulo, if I recall correctly.  
> **username_0**: Great! I've updated the python code to have all of the Ints where needed. I ran it again and it seems I'm still getting a negative number for the day value. But I'm sure it fixed some issues that were present before.  
> **username_0**: Just updated the code even further. I was accidentally printing out "f" when "dy" was supposed to be the printed as the first thing. Now I just need to fix the formatting so that I can use my "form" variable, and then double check output values.  
> **username_0**: The numbers still definitely aren't calculating right. I'm only getting one line of data for any year I enter. Also, the source is still not showing up.  

---

## 1 answers

### Answer by username_4 on 2015-11-24. Score: 2 | Accepted answer

There were a few other small errors in your translation (such as indenting too much stuff in two of the if statements), and I think you might have been trying to be too clever with changing SUBROUTINE into defined functions when it was a lot easier to just drop the code in-line since it was only used once.
There were minor tweaks and as far as I can tell this seems to work (you'll want to play with the formatted print to clean it up, and this runs under Python 3; if you need it to run under Python 2, you'll probably need to change the print and num1.write statements):
#!/usr/bin/env python

import math

print('Program to flag Jovian Decametric windows')
month = ['JAN','FEB','MAR','APR','MAY','JUN','JUL','AUG','SEP','OCT','NOV','DEC']
week = 42.46/360
pi = 3.141593
kr = pi / 180
num1 = open('jovrad.txt', 'w')
yy = int(input(("Year for which predictions are required ")))
e = math.trunc(((yy-1) / 100))
print(e)
f = 2 - e + math.trunc(e/4)
print(f)
jd = math.trunc(365.25 * (yy - 1)) + 1721423 + f + .5
print(jd)
d0 = jd - 2435108
print(d0)
incr = 0
dmax = 0
tx = 0
ty = 0
yyly = 0
yylc = 0
if yy / 400 - math.trunc((yy / 400)) == 0:
    incr = 1
    print("in if-1")
yyly = yy / 4 - math.trunc((yy / 4))
yylc = yy / 100 - math.trunc((yy / 100))
if yyly == 0 and yylc != 0:
    print("in if-2")
    incr = 1
ty = 59 + incr
dmax = 365 + incr
tx = ty + .5
num1.write("******************************************************\n")
pout = "  JOVIAN IO-DECAMETRIC EMISSION PREDICTIONS FOR " + str(yy) + "\n"
num1.write(pout)
num1.write("******************************************************\n")
num1.write("\n")
num1.write("Day   Date   Hr(UT)  Io_Phase   CML   Dist(AU)  Source")
num1.write("\n")
th = 0

while int(th / 24) + 1 <= dmax:
    d = d0 + (th / 24)
    v = (157.0456 + .0011159 * d) % 360
    m = (357.2148 + .9856003 * d) % 360
    n = (94.3455 + .0830853 * d + .33 * math.sin(kr * v)) % 360
    j = (351.4266 + .9025179 * d - .33 * math.sin(kr * v)) % 360
    a = 1.916 * math.sin(kr * m) + .02 * math.sin(kr * 2 * m)
    b = 5.552 * math.sin(kr * n) + .167 * math.sin(kr * 2 * n)
    k = j + a - b
    r = 1.00014 - .01672 * math.cos(kr * m) - .00014 * math.cos(kr * 2 * m)
    re = 5.20867 - .25192 * math.cos(kr * n) - .0061 * math.cos(kr * 2 * n)
    dt = math.sqrt(re * re + r * r - 2 * re * r * math.cos(kr * k))
    sp = r * math.sin(kr * k) / dt
    ps = sp / .017452
    dl = d - dt / 173
    pb = ps - b
    xi = 150.4529 * math.trunc((dl)) + 870.4529 * (dl - math.trunc((dl)))
    L3 = (274.319 + pb + xi + .01016 * 51) % 360
    U1 = 101.5265 + 203.405863 * dl + pb
    U2 = 67.81114 + 101.291632 * dl + pb
    z = (2 * (U1 - U2)) % 360
    U1 = U1 + .472 * math.sin(kr * z)
    U1 = (U1 + 180) % 360
    s = ""
    if L3 < 255 and L3 > 200 and U1 < 250 and U1 > 220:
        s = "Io-A"
    if L3 < 180 and L3 > 105 and U1 < 100 and U1 > 80:
        s = "Io-B"
    if L3 < 350 and L3 > 300 and U1 < 250 and U1 > 230:
        s = "Io-C"
    if s != "":
        dy = math.trunc((th / 24)) + 1
        h = th - (dy - 1) * 24
        if(dy > th):
            m = math.trunc((dy - tx) / 30.6) + 3
            da = dy - ty - math.trunc((m - 3) * 30.6 + .5)
        else:
            m = math.trunc((dy - 1) / 31) + 1
            da = dy - (m - 1) * 31
        mn = month[(m-1)]
#        mn = month[(m-1)*3+1:(m-1)*3-1+3]
        outstring = "%i  %s  %i  %5.3f  %5.3f  %5.3f  %5.3f  %s\n" % (dy, mn, da, h, U1, L3, dt, s)
        num1.write(outstring)
    th = th + .5

num1.close()
print("Program Complete  - results in file JOVRAD.TXT")

> **Comments:**
> **username_0**: Thank you so much! I really appreciate it. This does run 100% in Python 2 by the way. And yeah I'll have to work a little bit on the printing format.  
> **username_0**: Also, I removed the first variable in the "outstring" you wrote. It just starts with the month and then the day (variable "da").  

---

