/* Copyright (GPL) 2004 Mike Chirico mchirico@comcast.net
Updated: Sun Nov 28 15:15:05 EST 2004
Program adapted by Mike Chirico mchirico@comcast.net
Reference:
http://prdownloads.sourceforge.net/souptonuts/working_with_time.tar.gz?download
http://www.srrb.noaa.gov/highlights/sunrise/sunrise.html
Compile:
gcc -o sunrise -Wall -W -O2 -s -pipe -lm sunrise.c
or for debug output
gcc -o sunrise -DDEBUG=1 -Wall -W -O2 -s -pipe -lm sunrise.c
This can also go in a batch job to calculate the next
20 days as follows:
#!/bin/bash
lat=39.95
long=75.15
for (( i=0; i <= 20; i++))
do
./sunrise `date -d "+$i day" "+%Y %m %d"` $lat $long
done
*/
/* gcc -DDEBUG=1 .. */
#ifndef DEBUG
#define DEBUG 0
#endif
#include <sys/time.h>
#include <time.h>
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
double calcSunEqOfCenter(double t);
/* Convert degree angle to radians */
double degToRad(double angleDeg)
{
return (M_PI * angleDeg / 180.0);
}
double radToDeg(double angleRad)
{
return (180.0 * angleRad / M_PI);
}
double calcMeanObliquityOfEcliptic(double t)
{
double seconds = 21.448 - t*(46.8150 + t*(0.00059 - t*(0.001813)));
double e0 = 23.0 + (26.0 + (seconds/60.0))/60.0;
return e0; // in degrees
}
double calcGeomMeanLongSun(double t)
{
double L = 280.46646 + t * (36000.76983 + 0.0003032 * t);
while( (int) L > 360 )
{
L -= 360.0;
}
while( L < 0)
{
L += 360.0;
}
return L; // in degrees
}
double calcObliquityCorrection(double t)
{
double e0 = calcMeanObliquityOfEcliptic(t);
double omega = 125.04 - 1934.136 * t;
double e = e0 + 0.00256 * cos(degToRad(omega));
return e; // in degrees
}
double calcEccentricityEarthOrbit(double t)
{
double e = 0.016708634 - t * (0.000042037 + 0.0000001267 * t);
return e; // unitless
}
double calcGeomMeanAnomalySun(double t)
{
double M = 357.52911 + t * (35999.05029 - 0.0001537 * t);
return M; // in degrees
}
double calcEquationOfTime(double t)
{
double epsilon = calcObliquityCorrection(t);
double l0 = calcGeomMeanLongSun(t);
double e = calcEccentricityEarthOrbit(t);
double m = calcGeomMeanAnomalySun(t);
double y = tan(degToRad(epsilon)/2.0);
y *= y;
double sin2l0 = sin(2.0 * degToRad(l0));
double sinm = sin(degToRad(m));
double cos2l0 = cos(2.0 * degToRad(l0));
double sin4l0 = sin(4.0 * degToRad(l0));
double sin2m = sin(2.0 * degToRad(m));
double Etime = y * sin2l0 - 2.0 * e * sinm + 4.0 * e * y * sinm * cos2l0
- 0.5 * y * y * sin4l0 - 1.25 * e * e * sin2m;
return radToDeg(Etime)*4.0; // in minutes of time
}
double calcTimeJulianCent(double jd)
{
double T = ( jd - 2451545.0)/36525.0;
return T;
}
double calcSunTrueLong(double t)
{
double l0 = calcGeomMeanLongSun(t);
double c = calcSunEqOfCenter(t);
double O = l0 + c;
return O; // in degrees
}
double calcSunApparentLong(double t)
{
double o = calcSunTrueLong(t);
double omega = 125.04 - 1934.136 * t;
double lambda = o - 0.00569 - 0.00478 * sin(degToRad(omega));
return lambda; // in degrees
}
double calcSunDeclination(double t)
{
double e = calcObliquityCorrection(t);
double lambda = calcSunApparentLong(t);
double sint = sin(degToRad(e)) * sin(degToRad(lambda));
double theta = radToDeg(asin(sint));
return theta; // in degrees
}
double calcHourAngleSunrise(double lat, double solarDec)
{
double latRad = degToRad(lat);
double sdRad = degToRad(solarDec);
double HA = (acos(cos(degToRad(90.833))/(cos(latRad)*cos(sdRad))-tan(latRad) * tan(sdRad)));
return HA; // in radians
}
double calcHourAngleSunset(double lat, double solarDec)
{
double latRad = degToRad(lat);
double sdRad = degToRad(solarDec);
double HA = (acos(cos(degToRad(90.833))/(cos(latRad)*cos(sdRad))-tan(latRad) * tan(sdRad)));
return -HA; // in radians
}
double calcJD(int year,int month,int day)
{
if (month <= 2) {
year -= 1;
month += 12;
}
int A = floor(year/100);
int B = 2 - A + floor(A/4);
double JD = floor(365.25*(year + 4716)) + floor(30.6001*(month+1)) + day + B - 1524.5;
return JD;
}
double calcJDFromJulianCent(double t)
{
double JD = t * 36525.0 + 2451545.0;
return JD;
}
double calcSunEqOfCenter(double t)
{
double m = calcGeomMeanAnomalySun(t);
double mrad = degToRad(m);
double sinm = sin(mrad);
double sin2m = sin(mrad+mrad);
double sin3m = sin(mrad+mrad+mrad);
double C = sinm * (1.914602 - t * (0.004817 + 0.000014 * t)) + sin2m * (0.019993 - 0.000101 * t) + sin3m * 0.000289;
return C; // in degrees
}
double calcSunriseUTC(double JD, double latitude, double longitude)
{
double t = calcTimeJulianCent(JD);
// *** First pass to approximate sunrise
double eqTime = calcEquationOfTime(t);
double solarDec = calcSunDeclination(t);
double hourAngle = calcHourAngleSunrise(latitude, solarDec);
double delta = longitude - radToDeg(hourAngle);
double timeDiff = 4 * delta; // in minutes of time
double timeUTC = 720 + timeDiff - eqTime; // in minutes
double newt = calcTimeJulianCent(calcJDFromJulianCent(t) + timeUTC/1440.0);
eqTime = calcEquationOfTime(newt);
solarDec = calcSunDeclination(newt);
hourAngle = calcHourAngleSunrise(latitude, solarDec);
delta = longitude - radToDeg(hourAngle);
timeDiff = 4 * delta;
timeUTC = 720 + timeDiff - eqTime; // in minutes
return timeUTC;
}
double calcSunsetUTC(double JD, double latitude, double longitude)
{
double t = calcTimeJulianCent(JD);
// *** First pass to approximate sunset
double eqTime = calcEquationOfTime(t);
double solarDec = calcSunDeclination(t);
double hourAngle = calcHourAngleSunset(latitude, solarDec);
double delta = longitude - radToDeg(hourAngle);
double timeDiff = 4 * delta; // in minutes of time
double timeUTC = 720 + timeDiff - eqTime; // in minutes
double newt = calcTimeJulianCent(calcJDFromJulianCent(t) + timeUTC/1440.0);
eqTime = calcEquationOfTime(newt);
solarDec = calcSunDeclination(newt);
hourAngle = calcHourAngleSunset(latitude, solarDec);
delta = longitude - radToDeg(hourAngle);
timeDiff = 4 * delta;
timeUTC = 720 + timeDiff - eqTime; // in minutes
// printf("************ eqTime = %f \nsolarDec = %f \ntimeUTC = %f\n\n",eqTime,solarDec,timeUTC);
return timeUTC;
}
int main(int argc, char **argv)
{
time_t seconds;
time_t tseconds;
struct tm *ptm=NULL;
struct tm tm;
int year=2004,month=8,day=21,dst=-1;
char buffer[30];
float JD=calcJD(year,month,day);
double latitude = 39.95; //convert to just degrees. No min/sec
double longitude = 75.15;
if(argc == 6) {
year=atoi(argv[1]);
month=atoi(argv[2]);
day=atoi(argv[3]);
latitude = atof(argv[4]);
longitude = atof(argv[5]);
JD = calcJD(year,month,day);
} else {
printf("\nUsage: (note convert latitude and longitude to degrees) \n");
printf("./sunrise year month date latitude longitude \n\n");
printf(" Latitude and longitude must be in degrees and or fraction of degrees. Not min sec\n");
printf("\n US Listings Of Latitude and Longitude\n");
printf(" http://www.census.gov/geo/www/tiger/latlng.txt\n");
printf(" But use positive values for longitude for above listing\n");
printf("\n Example: (Just outside Philadelphia PA, USA)\n\n");
printf("./sunrise 2004 8 21 39.95 75.15 \n");
printf("./sunrise `date \"+%cY %cm %cd\"` 39.95 75.15\n\n",'%','%','%');
}
if(DEBUG)
printf("Julian Date %f \n",JD);
if(DEBUG)
printf("Sunrise timeUTC %lf \n", calcSunriseUTC( JD, latitude, longitude));
if(DEBUG)
printf("Sunset timeUTC %lf \n", calcSunsetUTC( JD, latitude, longitude));
tm.tm_year= year-1900;
tm.tm_mon=month-1; /* Jan = 0, Feb = 1,.. Dec = 11 */
tm.tm_mday=day;
tm.tm_hour=0;
tm.tm_min=0;
tm.tm_sec=0;
tm.tm_isdst=-1;
seconds = mktime(&tm);
if(DEBUG)
printf("Number of seconds %ld\n",seconds);
if(DEBUG)
printf("Number of year %d\n",year);
int delta;
dst=tm.tm_isdst;
ptm = gmtime ( &seconds );
delta= ptm->tm_hour;
if(DEBUG)
printf("delta=%d dst=%d\n",delta,dst);
tseconds= seconds;
if(DEBUG)
printf("Number of seconds %ld\n",seconds);
seconds= seconds + calcSunriseUTC( JD, latitude, longitude)*60;
seconds= seconds - delta*3600;
strftime(buffer,30,"%m-%d-%Y %T",localtime(&seconds));
printf("Sunrise %s ",buffer);
seconds=tseconds;
seconds+=calcSunsetUTC( JD, latitude, longitude)*60;
seconds= seconds - delta*3600;
strftime(buffer,30,"%m-%d-%Y %T",localtime(&seconds));
printf("Sunset %s\n",buffer);
return 0;
}
Linux System Admin Tips: There are over 200 Linux tips and tricks in this article. That is over 100 pages covering everything from NTP, setting up 2 IP address on one NIC, sharing directories among several users, putting running jobs in the background, find out who is doing what on your system by examining open sockets and the ps command, how to watch a file, how to prevent even root from deleting a file, tape commands, setting up cron jobs, using rsync, using screen conveniently with emacs, how to kill every process for a user, security tips and a lot more. These tip grow weekly. The above link will download the text version for easy grep searching. There is also an html version here.
Breaking Firewalls with OpenSSH and PuTTY: If the system administrator deliberately filters out all traffic except port 22 (ssh), to a single server, it is very likely that you can still gain access other computers behind the firewall. This article shows how remote Linux and Windows users can gain access to firewalled samba, mail, and http servers. In essence, it shows how openSSH and Putty can be used as a VPN solution for your home or workplace.
MySQL Tips and Tricks: Find out who is doing what in MySQL and how to kill the process, create binary log files, connect, create and select with Perl and Java, remove duplicates in a table with the index command, rollback and how to apply, merging several tables into one, updating foreign keys, monitor port 3306 with the tcpdump command, creating a C API, complex selects, and much more.
Create a Live Linux CD - BusyBox and OpenSSH Included: These steps will show you how to create a functioning Linux system, with the latest 2.6 kernel compiled from source, and how to integrate the BusyBox utilities including the installation of DHCP. Plus, how to compile in the OpenSSH package on this CD based system. On system boot-up a filesystem will be created and the contents from the CD will be uncompressed and completely loaded into RAM -- the CD could be removed at this point for boot-up on a second computer. The remaining functioning system will have full ssh capabilities. You can take over any PC assuming, of course, you have configured the kernel with the appropriate drivers and the PC can boot from a CD. This tutorial steps you through the whole processes.
SQLite Tutorial : This article explores the power and simplicity of sqlite3, first by starting with common commands and triggers, then the attach statement with the union operation is introduced in a way that allows multiple tables, in separate databases, to be combined as one virtual table, without the overhead of copying or moving data. Next, the simple sign function and the amazingly powerful trick of using this function in SQL select statements to solve complex queries with a single pass through the data is demonstrated, after making a brief mathematical case for how the sign function defines the absolute value and IF conditions.
The Lemon Parser Tutorial: This article explains how to build grammars and programs using the lemon parser, which is faster than yacc. And, unlike yacc, it is thread safe.
How to Compile the 2.6 kernel for Red Hat 9 and 8.0 and get Fedora Updates: This is a step by step tutorial on how to compile the 2.6 kernel from source.
Virtual Filesystem: Building A Linux Filesystem From An Ordinary File. You can take a disk file, format it as ext2, ext3, or reiser filesystem and then mount it, just like a physical drive. Yes, it then possible to read and write files to this newly mounted device. You can also copy the complete filesystem, since it is just a file, to another computer. If security is an issue, read on. This article will show you how to encrypt the filesystem, and mount it with ACL (Access Control Lists), which give you rights beyond the traditional read (r) write (w) and execute (x) for the 3 user groups file, owner and other.
Working With Time: What? There are 61 seconds in a minute? We can go back in time? We still tell time by the sun?
Mike Chirico, a father of triplets (all girls) lives outside of
Philadelphia, PA, USA. He has worked with Linux since 1996, has a Masters
in Computer Science and Mathematics from Villanova University, and has
worked in computer-related jobs from Wall Street to the University of
Pennsylvania. His hero is Paul Erdos, a brilliant number theorist who was
known for his open collaboration with others.
Mike's notes page is souptonuts. For
open source consulting needs, please send an email to
mchirico@gmail.com. All consulting work must include a donation to
SourceForge.net.
