I understood generally what the code was doing from the variable names, but I don't know why you wanted to do those things. When writing software, I believe that comments should describe why you're doing something rather than describe what the code is doing because the code itself always says what it's doing.

I'm interested in the meanings of the physical quantities behind the terms "solar declination" and "angle of incidence" as used here and why you are using them in the manner shown. When I did these sorts of calculations last century, I used some simplifying assumptions like a circular orbit, constant orbital phase, and a day much shorter than the orbital period to get some very simple results. I came up with a spreadsheet that gave the set of graphs below for the traditional seasons:

Click image for larger version.  Name: insolation.png  Views: 10  Size: 53.5 KB  ID: 135196

They aren't "accurate" by any stretch of the imagination, but certainly were good enough for my purposes at the time, especially considering the coarse level of fidelity in what I was doing elsewhere. Using Earth as an example, the 4% variance in energy between perihelion and aphelion is dwarfed by the error of sampling the orbit at only four points.

As far as tooling, I used a spreadsheet for this because the calculations are mostly a function of two values (tilt and latitude), which is really easy to do on a spreadsheet. And it's what I had at hand way back when, of course.