Satellite Observers Workbench. NOT yet complete, just published for forum posters to \"cherry pick\" pieces of code as requiered as an example.
satapi.c
00001 /**************************************************************************** 00002 * Copyright 2010 Andy Kirkham, Stellar Technologies Ltd 00003 * 00004 * This file is part of the Satellite Observers Workbench (SOWB). 00005 * 00006 * SOWB is free software: you can redistribute it and/or modify 00007 * it under the terms of the GNU General Public License as published by 00008 * the Free Software Foundation, either version 3 of the License, or 00009 * (at your option) any later version. 00010 * 00011 * SOWB is distributed in the hope that it will be useful, 00012 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00013 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00014 * GNU General Public License for more details. 00015 * 00016 * You should have received a copy of the GNU General Public License 00017 * along with SOWB. If not, see <http://www.gnu.org/licenses/>. 00018 * 00019 * $Id: main.cpp 5 2010-07-12 20:51:11Z ajk $ 00020 * 00021 ***************************************************************************/ 00022 00023 #include "sowb.h" 00024 #include "user.h" 00025 #include "satapi.h" 00026 #include "utils.h" 00027 #include "debug.h" 00028 #include "gpio.h" 00029 #include "osd.h" 00030 #include "nexstar.h" 00031 #include "utils.h" 00032 00033 #ifndef M_PI 00034 #define M_PI 3.1415926535898 00035 #endif 00036 00037 #define SCAN_INTERVAL 60 00038 00039 double satapi_aos(SAT_POS_DATA *q, bool goto_aos) { 00040 double tsince; 00041 char temp1[32], temp2[32]; 00042 00043 strcpy(q->elements[0], "Lacrosse 2"); 00044 strcpy(q->elements[1], "1 21147U 91017A 10269.82093092 0.00000020 00000-0 28786-5 0 05"); 00045 strcpy(q->elements[2], "2 21147 67.9820 220.2995 0002000 244.7811 115.2189 14.76261286 07"); 00046 00047 observer_now(q); 00048 00049 P22_ASSERT; 00050 00051 for (q->tsince = 0; q->tsince < (SCAN_INTERVAL * 90); q->tsince += SCAN_INTERVAL) { 00052 KICK_WATCHDOG; /* We are busy! */ 00053 satallite_calculate(q); 00054 if (q->elevation >= 10.) { 00055 /* Above horizon viewing. Work back to AOS. */ 00056 for (tsince = q->tsince, q->tsince--; q->elevation > 10. ; tsince = q->tsince--) { 00057 satallite_calculate(q); 00058 if (q->elevation < 10.) { 00059 //sprintf(temp, "%03f Q AOS El:%.1f AZ:%.1f %dKm\r\n", q->tsince, q->elevation, q->azimuth, (int)q->range); 00060 //debug_printf(temp); 00061 q->tsince = tsince; 00062 satallite_calculate(q); 00063 sprintf(temp1, "%03f T AOS El:%.1f AZ:%.1f %dKm\r\n", q->tsince, q->elevation, q->azimuth, (int)q->range); 00064 debug_printf(temp1); 00065 P22_DEASSERT; 00066 if (goto_aos) { 00067 sprintf(temp1, "%s T-%.2f", q->elements[0], tsince); 00068 osd_string_xy(1, 12, temp1); 00069 sprintf(temp1, "AOS %.2f%c %s%c %dKm", q->elevation, 176, printDouble_3_2(temp2, q->azimuth), 176, (int)q->range); 00070 osd_string_xy(1, 13, temp1); 00071 _nexstar_goto((uint32_t)((q->elevation / 360.) * 65536), (uint32_t)((q->azimuth / 360.) * 65536)); 00072 } 00073 return tsince; 00074 } 00075 } 00076 } 00077 } 00078 00079 P22_DEASSERT; 00080 return 0.; 00081 } 00082 00083 int satallite_calculate(SAT_POS_DATA *q) { 00084 double tsince; 00085 00086 /* Ensure the time and place are valid. */ 00087 if (!q->time.is_valid) return -1; 00088 if (!q->location.is_valid) return -2; 00089 00090 ClearFlag(ALL_FLAGS); 00091 00092 Get_Next_Tle_Set(q->elements, &q->tle); 00093 00094 select_ephemeris(&q->tle); 00095 00096 q->jd_utc = gps_julian_date(&q->time); 00097 q->jd_epoch = Julian_Date_of_Epoch(q->tle.epoch); 00098 00099 tsince = ((q->jd_utc + (q->tsince * (1 / 86400.))) - q->jd_epoch) * xmnpda; 00100 00101 if (isFlagSet(DEEP_SPACE_EPHEM_FLAG)) { 00102 SDP4(tsince, &q->tle, &q->pos, &q->vel, &q->phase); 00103 } 00104 else { 00105 SGP4(tsince, &q->tle, &q->pos, &q->vel, &q->phase); 00106 } 00107 00108 Convert_Sat_State(&q->pos, &q->vel); 00109 SgpMagnitude(&q->vel); // scalar magnitude, not brightness... 00110 q->velocity = q->vel.w; 00111 00112 /* Populate the geodetic_t struct from data supplied. */ 00113 q->observer.lat = q->location.latitude * de2ra; 00114 q->observer.lon = q->location.longitude * de2ra; 00115 q->observer.alt = q->location.height / 1000.; 00116 if (q->location.north_south == 'S') q->observer.lat *= -1.; 00117 if (q->location.east_west == 'W') q->observer.lon *= -1.; 00118 00119 Calculate_Obs(q->jd_utc, &q->pos, &q->vel, &q->observer, &q->obs_set); 00120 Calculate_LatLonAlt(q->jd_utc, &q->pos, &q->sat_geodetic); 00121 00122 q->azimuth = Degrees(q->obs_set.x); 00123 q->elevation = Degrees(q->obs_set.y); 00124 q->range = q->obs_set.z; 00125 q->rangeRate = q->obs_set.w; 00126 q->height = q->sat_geodetic.alt; 00127 00128 return 0; 00129 } 00130 00131 /** observer_now 00132 * 00133 * Fills the data structure with the observers time and position. 00134 * 00135 * @param SAT_POS_DATA * A pointer to the data structure. 00136 */ 00137 SAT_POS_DATA * observer_now(SAT_POS_DATA *q) { 00138 gps_get_time(&q->time); 00139 gps_get_location_average(&q->location); 00140 return q; 00141 } 00142 00143 AltAz * radec2altaz(double siderealDegrees, GPS_LOCATION_AVERAGE *location, RaDec *radec, AltAz *altaz) { 00144 double HA, DEC, LAT, mul, altitude, azimuth; 00145 00146 mul = location->north_south == 'S' ? -1.0 : 1.0; 00147 00148 /* Convert to radians. */ 00149 HA = siderealDegrees * (M_PI / 180.0) - (radec->ra * (M_PI / 180)); 00150 DEC = radec->dec * (M_PI / 180.0); 00151 LAT = (location->latitude * mul) * (M_PI / 180.0); 00152 00153 altitude = atan2(- sin(HA) * cos(DEC), cos(LAT) * sin(DEC) - sin(LAT) * cos(DEC) * cos(HA)); 00154 azimuth = asin(sin(LAT) * sin(DEC) + cos(LAT) * cos(DEC) * cos(HA)); 00155 00156 // Convert to degrees and swing azimuth around if needed. 00157 altaz->alt = azimuth * 180.0 / M_PI; 00158 altaz->azm = altitude * 180.0 / M_PI; 00159 if (altaz->azm < 0) altaz->azm += 360.0; 00160 00161 return altaz; 00162 } 00163 00164 RaDec * altaz2radec(double siderealDegrees, GPS_LOCATION_AVERAGE *location, AltAz *altaz, RaDec *radec) { 00165 double ALT, AZM, LAT, HA, DEC, mul; 00166 00167 mul = location->north_south == 'S' ? -1.0 : 1.0; 00168 00169 /* Convert to radians. */ 00170 LAT = (location->latitude * mul) * (M_PI / 180.0); 00171 ALT = altaz->alt * (M_PI / 180.0); 00172 AZM = altaz->azm * (M_PI / 180.0); 00173 00174 /* Calculate the declination. */ 00175 DEC = asin( ( sin(ALT) * sin(LAT) ) + ( cos(ALT) * cos(LAT) * cos(AZM) ) ); 00176 radec->dec = DEC * 180.0 / M_PI; 00177 while (radec->dec < 0.0) radec->dec += 360.0; 00178 while (radec->dec > 360.0) radec->dec -= 360.0; 00179 00180 /* Calculate the hour angle. */ 00181 HA = ( acos((sin(ALT) - sin(LAT) * sin(DEC)) / (cos(LAT) * cos(DEC)))) * 180.0 / M_PI; 00182 if (sin(AZM) > 0.0) HA = 360.0 - HA; 00183 00184 /* Correct the HA for our sidereal time. */ 00185 HA = (siderealDegrees / 360.0 * 24.0) - (HA / 15.0); 00186 if (HA < 0.0) HA += 24.0; 00187 00188 /* Convert the HA into degrees for the return. */ 00189 radec->ra = HA / 24.0 * 360.0; 00190 00191 return radec; 00192 } 00193
Generated on Tue Jul 12 2022 18:05:35 by 1.7.2