1
0
Fork 0
pull/2/head
Cees Bassa 2013-09-20 16:57:15 +02:00
parent bdb395e139
commit 6098ee7d86
4 changed files with 142 additions and 14 deletions

View File

@ -68,11 +68,11 @@ int main(int argc,char *argv[])
img=read_fits(argv[1],0);
// Hard coded
img.ra0=209;
img.de0=18.00;
sx=-30.0;
sy=30.0;
q=-60.0;
img.ra0=323.5;
img.de0=6.19;
sx=-10.0;
sy=10.0;
q=-00.0;
img.x0=0.5*(float) img.naxis1;
img.y0=0.5*(float) img.naxis2;

View File

@ -13,7 +13,9 @@ struct image {
float xmin,xmax,ymin,ymax;
double tr[6];
double mjd;
char nfd[32];
char nfd[32],observer[32];
int cospar;
float exptime;
};
struct image read_jpg(char *filename);
void write_fits(struct image img,char *filename);
@ -27,9 +29,12 @@ int main(int argc,char *argv[])
struct image img;
char infile[64],outfile[64],nfd[32];
double mjd,delay=0.0,tz=0.0;
int cospar=0;
char observer[32]="Cees Bassa";
float exptime=10.06;
// Decode options
while ((arg=getopt(argc,argv,"i:t:o:d:Z:"))!=-1) {
while ((arg=getopt(argc,argv,"i:t:o:d:Z:c:T:O:"))!=-1) {
switch(arg) {
case 'i':
@ -52,6 +57,18 @@ int main(int argc,char *argv[])
strcpy(nfd,optarg);
break;
case 'c':
cospar=atoi(optarg);
break;
case 'T':
exptime=atof(optarg);
break;
case 'O':
strcpy(observer,optarg);
break;
default:
return 0;
@ -72,6 +89,11 @@ int main(int argc,char *argv[])
img.mjd=mjd;
}
// Set properties
img.cospar=cospar;
img.exptime=exptime;
strcpy(img.observer,observer);
if (outfile!=NULL)
write_fits(img,outfile);
@ -227,9 +249,12 @@ void write_fits(struct image img,char *filename)
qfits_header_add(qh,"DATE-OBS",val," ",NULL);
sprintf(val,"%lf",img.mjd);
qfits_header_add(qh,"MJD-OBS",val," ",NULL);
qfits_header_add(qh,"COSPAR","0000"," ",NULL);
qfits_header_add(qh,"EXPTIME","5.00"," ",NULL);
qfits_header_add(qh,"OBSERVER","Cees Bassa"," ",NULL);
sprintf(val,"%d",img.cospar);
qfits_header_add(qh,"COSPAR",val," ",NULL);
sprintf(val,"%f",img.exptime);
qfits_header_add(qh,"EXPTIME",val," ",NULL);
sprintf(val,"%s",img.observer);
qfits_header_add(qh,"OBSERVER",val," ",NULL);
// Dump fitsheader
// qfits_header_dump(qh,stdout);

View File

@ -51,6 +51,28 @@ struct image read_fits(char *filename,int pnum);
int fgetline(FILE *file,char *s,int lim);
int select_nearest(struct catalog c,float x,float y);
// Return x modulo y [0,y)
double modulo(double x,double y)
{
x=fmod(x,y);
if (x<0.0) x+=y;
return x;
}
// Greenwich Mean Sidereal Time
double gmst(double mjd)
{
double t,gmst;
t=(mjd-51544.5)/36525.0;
gmst=modulo(280.46061837+360.98564736629*(mjd-51544.5)+t*t*(0.000387933-t/38710000),360.0);
return gmst;
}
void plot_objects(char *filename)
{
int i;
@ -336,7 +358,7 @@ void reduce_point(struct observation *obs,struct image img,float tmid,float x,fl
int iframe,k;
double ra,de,rx,ry;
float dx,dy,dt;
double mjd;
double mjd,mjd1,mjd2;
char nfd[32],sra[15],sde[15];
// Transform position
@ -346,6 +368,11 @@ void reduce_point(struct observation *obs,struct image img,float tmid,float x,fl
ry=img.b[0]+img.b[1]*dx+img.b[2]*dy;
reverse(img.ra0,img.de0,rx,ry,&ra,&de);
// Correct for stationary camera
mjd1=img.mjd+0.5*(double) img.exptime/86400.0;
mjd2=img.mjd+(double) tmid/86400.0;
ra+=gmst(mjd2)-gmst(mjd1);
dec2sex(ra/15.0,sra,0);
dec2sex(de,sde,1);

View File

@ -30,6 +30,7 @@ struct map {
char orientation[LIM],projection[4],observer[32];
char nfd[LIM],starfile[LIM],tlefile[LIM],iodfile[LIM],xyzfile[LIM];
char datadir[LIM],tledir[LIM];
float saltmin;
double lat,lng;
double h,sra,sde,sazi,salt;
float alt,timezone;
@ -102,12 +103,15 @@ void get_site(int site_id);
void usage()
{
printf("skymap t:c:i:R:D:hs:d:l:P:r:V:p:\n\n");
printf("skymap t:c:i:R:D:hs:d:l:P:r:V:p:A:E:S:\n\n");
printf("t date/time (yyyy-mm-ddThh:mm:ss.sss) [default: now]\n");
printf("c TLE catalog file [default: classfd.tle]\n");
printf("i satellite ID (NORAD) [default: all]\n");
printf("R R.A.\n");
printf("D Decl.\n");
printf("A Azimuth\n");
printf("E Elecation\n");
printf("S All night\n");
printf("h this help\n");
printf("s site (COSPAR)\n");
printf("d IOD observations\n");
@ -153,6 +157,7 @@ void init_skymap(void)
m.visflag=0;
m.planar=0;
m.agelimit=-1.0;
m.saltmin=-6.0;
// Default settings
strcpy(m.observer,"Unknown");
@ -405,6 +410,59 @@ void plot_xyz(double mjd0,char *filename)
return;
}
void allnight(void)
{
int flag;
xyz_t sunpos;
double ra,de,azi,alt,alt0;
double mjd,mjdrise=-1.0,mjdset=-1.0;
char nfd[32];
// Find first event
for (flag=0,mjd=m.mjd;mjd<m.mjd+1.0;mjd+=1.0/86400) {
sunpos_xyz(mjd,&sunpos,&ra,&de);
equatorial2horizontal(mjd,ra,de,&azi,&alt);
if (flag!=0) {
if (alt>m.saltmin && alt0<=m.saltmin)
mjdrise=mjd;
if (alt<m.saltmin && alt0>=m.saltmin)
mjdset=mjd;
}
if (flag==0)
flag=1;
alt0=alt;
}
if (mjdrise<mjdset) {
for (flag=0,mjd=mjdrise-1.0;mjd<mjdrise;mjd+=1.0/86400) {
sunpos_xyz(mjd,&sunpos,&ra,&de);
equatorial2horizontal(mjd,ra,de,&azi,&alt);
if (flag!=0) {
if (alt>m.saltmin && alt0<=m.saltmin)
mjdrise=mjd;
if (alt<m.saltmin && alt0>=m.saltmin)
mjdset=mjd;
}
if (flag==0)
flag=1;
alt0=alt;
}
}
m.mjd=mjdset;
mjd2date(m.mjd,m.nfd);
mjd2date(mjdrise,nfd);
m.length=(mjdrise-mjdset)*86400;
printf("Set: %s; Rise: %s; length: %.0fs\n",m.nfd,nfd,m.length);
return;
}
int main(int argc,char *argv[])
{
int i,arg=0;
@ -415,7 +473,7 @@ int main(int argc,char *argv[])
init_skymap();
// Decode options
while ((arg=getopt(argc,argv,"t:c:i:R:D:hs:d:l:P:r:V:p:"))!=-1) {
while ((arg=getopt(argc,argv,"t:c:i:R:D:hs:d:l:P:r:V:p:A:E:S:"))!=-1) {
switch(arg) {
case 't':
@ -424,6 +482,11 @@ int main(int argc,char *argv[])
m.iodpoint=-1;
break;
case 'S':
m.saltmin=atof(optarg);
allnight();
break;
case 'c':
strcpy(m.tlefile,optarg);
break;
@ -485,6 +548,19 @@ int main(int argc,char *argv[])
m.level=5;
break;
case 'A':
m.azi0=modulo(atof(optarg)+180.0,360.0);
strcpy(m.orientation,"horizontal");
m.level=3;
break;
case 'E':
m.alt0=atof(optarg);
strcpy(m.orientation,"horizontal");
m.level=3;
break;
case 'h':
usage();
return 0;
@ -497,7 +573,7 @@ int main(int argc,char *argv[])
}
init_plot("/xs",10,0.75);
plot_skymap();
cpgend();