Sunday, July 12, 2009

C++ question?

How can I change the type of double in the vector class to the struct track from the following code


struct Track


{


string trackName;





int trackDuration;


};





class Vector


{


private:


double *array;


int num_elements;


int max;





// doubles the capacity of the vector by allocating a new array and copying


// all elements over to it


// used by the push_back function when the vector becomes full


void grow();





public:


// Constructors to initialize the vector's member variables





Vector() { max = 10; array = new double[max]; num_elements = 0; }


// default constructor





Vector(int, double);


// initialize all elements to a specific value





// destructor


~Vector() { delete[] array; }





// Copy Constructor -- Added 4/30/07


Vector(const Vector%26amp;);





// Assignment operator -- Added 4/30/07


Vector%26amp; operator=(const Vector%26amp;);





void push_back(double); //

C++ question?
Astronomy C/C++ source code





Last updated 5 January 2007





From time to time, I'm asked to provide source code for doing some sort of astronomical calculation, or for providing direct access to the numerous compressed datasets on the Guide CD-ROM. In the future, I intend to gather information about and links to all such source code on this page.





If you don't see what you want here, please ask. I've got a fair bit of source code not yet posted simply because I've not gotten around to documenting it thoroughly yet.





The source code in question is written in C; click here for info on source code in other languages, or here for info on other C/C++ source code for astronomy. I use an extension of .CPP just to force the compiler to do certain error checking (type casts, etc.) that C++ compilers do better than C compilers. The only use I make of C++ is in Windows user interface (MFC) code, and I haven't posted any of that yet. However, Mark Huss has produced a C++ implementation of much of this code, and you can click here to download the ZIPped C++ version (about 128 KBytes). If that link fails, try this one.





I will post remaining code here from time to time. If you have any priorities (bits you would really like to see soon), please let me know; I can probably rearrange the order to get your code uploaded faster.





Changes: Assorted changes have been made in the past with the regrettable absence of a change log. In the future, the .ZIP file will contain a CHANGES.TXT file, explaining what's different about the newly-posted files.





Porting: As much as possible, I've stuck with ANSI C. Most of the code has already been used in 32-bit DOS and Windows; some has even run in 16-bit DOS and Windows. A make file for Linux is provided, and the code works there (but I've only used it to create a Linux version of Find_Orb, and not much else.) Mark Huss' C++ version also compiles in Linux. The only warning I'd offer in this regard is that I haven't dealt with a "wrong-end" byte order platform yet, and I am quite certain that some of the code will crash in grisly fashion on such systems.





Documentation: Almost all the code has quite thorough descriptions of how the innards work. (This requirement is the main one slowing posting of code. Some source has been posted for Find_Orb and Charon. But most of it is complex enough that posting it without explanation would not be very helpful.)





Copyright restrictions: All this code is free for use in non-commercial applications. If you wish to use any of it in a commercial application, please let me know; I'm not averse to a little bartering in such cases.





There are two exceptions to this: the RealSky/DSS image extraction software and the FITS image compression software are entirely in the public domain. (Certain parts of the RealSky/DSS code came from the Space Telescope Science Institute (STScI), and the FITS compression software is a very slightly modified version of an STScI package.)





For all the source code, I would ask that you inform me of any bugs you find, and if you make interesting improvements, I'd very much like to hear about them.





Astronomy source code in other languages (so far, only Java): Personally, I haven't written much of anything in a language other than C or C++ since 1988. (Not a good thing, I realize.) I get occasional questions about source code in Java and Pascal (surprisingly, not in BASIC... though you can click here for a list of BASIC programs from Sky %26amp; Telescope.)





In addition to his C++ astronomy source code, Mark Huss has written a lot of astronomy source code in Java, and has updated and improved it a bit recently (March 2002) as described in the following e-mail:





hi Bill,





I finally got things organized - the updated library, source, and info is


now available on my website: http://mhuss.com/AstroLib.html.





javadoc is an amazing thing - take a glance at


http://mhuss.com/AstroLib/docs/api/index...





regards,





--mark huss





(And from a previous e-mail...)





This is what's currently running the lunar phase info and 'darkest hours'


page on the DVAA website (http://dvaa.org). You should know, however, that I


used Meeus 'mystic formulas' instead of applying your 'witness the truth


thereof' philosophy in a couple of places for expediency ;-).





regards,





--mark





Some other implementation of assorted astronomical algorithms in C/C++: Mark Huss has written quite a bit of astronomy-related code in C++. P. J. Naughter has implemented much of Meeus' Astronomical Algorithms (second edition) in C++. You can click here for his code and info about it.





Source code for accessing Guide datasets and orbital elements: There is rather a lot of code in this category, and it has its own separate page. Click here for information on this subject. There are also some pieces of code to access some large catalogues here.





Source code for accessing JPL's DE ephemerides: This has grown to the point where it needed a page of its own. Click here for C/C++ source code for accessing JPL DE ephemerides.





Source code for RealSky/DSS image extraction: This is the same library used in Guide for image extraction; it's also used in version of SkyMap Pro, and, I suspect, in some other commercial astronomy software. You can read about it and download it from this page.





Source code for FITS image compression: A very slightly modified version of an STScI program. Click here for information and source code.





Source code and a DOS executable for tables of Galilean satellite events: This is the code I used to generate the lists of occultations, eclipses, transits, and shadow events shown in Guide. Click here for an example of the output of this program. The code is quite ugly, and I don't recommend its use very highly. Click here to download the source code and .EXE files (about 69 KBytes).





Source code for basic astronomical calculations: This file will eventually encompass all the low-level routines I use in Guide: nutation, aberration, computing planet positions from VSOP and PS_1996, refraction, lunar positions, the orientations of planets, coordinate systems, sidereal time...





At present, it has many of these, but not all. Here's a summary of the purpose of most of the .CPP files:





Click here to download the ZIPped code (about 93 KBytes). Here's some discussion of the meaning of the various source files supplied:





ALT_AZ.CPP: Code to convert RA/dec to alt/az, and galactic coordinates to RA/dec.





ASTEPHEM.CPP: The main function for the code to compute asteroid ephemerides. It demonstrates accessing orbital elements from the Guide disk, and then shows how to compute an RA/dec, magnitude, elongation, and other data based on using those elements to compute the asteroid position, and VSOP to compute the Earth's position. The most significant thing omitted right now is the topocentric offset.





ASTFUNCS.CPP: Some basic functions to set up asteroid/comet elements and to compute their locations (solving Kepler's equation and so forth.)





The Kepler's equation solution draws quite a bit of interest, so much so that there is a separate page discussing Kepler's Equation.





BIG_VSOP.CPP: Code to compute planetary positions from the full VSOP series, which is stored in binary form on the Guide CD-ROM. The logic is very similar to that used for the "truncated" VSOP (see VSOPSON.CPP), except that the function works by reading in pieces of a file; VSOPSON.CPP assumes the data file has been read into memory.





CLASSEL.CPP: Takes a state vector (position and velocity) and computes the "classical elements" (semimajor axis, eccentricity, mean anomaly, argument of perihelion, etc.) Currently, the FIND_ORB orbit determination software and the integrat numerical integration program make use of this; after computing an orbit in vector format, they use this routine to produce something in a form suitable for human consumption.





COLORS.CPP and COLORS2.CPP: Given, say, a B-V color value, it's possible to compute a "pretty good" V-R value. Similar conversions are possible between most other color indices, with varying degrees of accuracy. The general form is a polynomial, say,





index2 = a0 + a1 * index + a2 * index2 + a3 * index3 +....


These resemble Taylor series, though I suspect they are actually curve-fit polynomials. I got them from some Fortran snippets supplied by Brian Skiff, and they have been converted to C here.





Certain conversions are handled by inverting the polynomial. The need to do this indicates a very unreliable color conversion, to be used only in cases of dire emergency.





Also, code to convert Tycho B and V to Johnson B and V is given, along with a little piece of test code. (That Tycho-to-Johnson algorithm is one based the book that was distributed with the Hipparcos data, "Introduction and Guide to the Data." COLORS2.CPP is a more recent, and probably better, Tycho-to-Johnson piece of source code.)





COM_FILE.CPP: This is pretty much Guide-specific code. Guide goes through some odd hoops to draw comets as speedily as possible. One headache it has to deal with is that there isn't any consolidated set of comet orbital elements, the way there is with asteroids (e.g., MPCORB or the Lowell Observatory ASTORB). Instead, it has to dynamically combine "current comet" data from MPC with a list of historical comet data. The confused nature of comet designations makes that a pain (it's not always easy to match comets from one list with those from another). Also, Guide uses some pregenerated data to omit comets that are too faint to be of interest over a desired time span... all nifty stuff, and it helps Guide to run somewhat more briskly. But it's probably not of interest to most people.





COSPAR.CPP: Code to compute the "planet orientation" matrix for the Sun, all planets, and quite a few satellites. A few years ago, the IAU set up a committee to define pole positions and longitude systems for these objects. The pole position is usually defined to be a linear function in RA/dec, with the meridian intersecting the J2000 plane as another linear function. Often, some small trig terms are added to those linear functions.





At present, COSPAR.CPP includes many, but not all, of the satellites shown by Guide. Some (such as the "captured asteroids" of gas giants) have no defined orientations. Others (small inner objects such as Amalthea) just haven't been tackled by me yet (no big job to do, but not much interest, either).





DATE.CPP: Source code for the following calendrical systems: Julian, Gregorian, Hebrew, Persian (Jalaali and 'Modern Persian'), Islamic, and Chinese, and French Republican. Basically, you can do conversions between JD and the calendar systems.





DELTA_T.CPP: Code to compute Delta_T = TD - UT for any date, using a lookup table for years 1620 to 2002 and extrapolations beyond those ranges.





DE_PLAN.CPP: Code to compute very precise planetary positions using the PS-1996 series of Chapront. Includes Mercury through Pluto. To use this, you'll need to download this file of coefficients.





DIST_PA.CPP: Code to compute the distance and position angle between two spherical coordinates. Theoretically, this ought to be a straightforward task. In practice, there are an amazing number of ways to get roundoff/truncation errors and divides by zero, mostly involving how you deal with cases where the distance between the two points is small. I am reasonably sure now that this code evades all such cases.





EART2000.CPP: Code to compute the Earth's location relative to the Sun, in J2000 coordinates, using VSOP.





EASTER.CPP: Code to figure out the date of Easter for a given (Gregorian) year.





ELP82DAT.CPP: Code to compute the position of the moon using the ELP-82 (Ephemerides Lunaire Parisienne 1982) theory of Chapront and Chapront-Touzé. To make use of this code, you'll need to download this file of coefficients.





GETPLANE.CPP: Code to compute the location of any planet, or the Moon, relative to the Sun. The resulting vector is given in five different systems.





JD.CPP: Produces an "example" program showing the use of the functions in DATE.CPP. You can type in, say, "jd 10 7 1999", and be told what day 10 July 1999 corresponds to in assorted calendars. (You can also get the actual JD corresponding to that date, or type in, say, "JD 2451324.879" and get the calendar date... which was the original reason I wrote the program.)





JSATS.CPP: Code to compute the positions of the Galilean satellites, using Lieske's theory as given in Meeus' Astronomical Algorithms. Quite a few smaller terms are omitted.





LUNAR2.CPP: Code to compute the position of the moon, using the truncated series from Meeus' Astronomical Algorithms. This uses the same data file (VSOP.BIN) as the VSOPSON.CPP code discussed below.





MISCELL.CPP: Code for basic matrix functions (inverting an orthogonal matrix, rotation, setting identity matrices, etc.) Also code to handle the odd system used for variable star designations, and to provide a version of the C-language ctime() function that takes a JD and can work over a full time range (ctime only works between 1970 and 2028.) Also, this version gives plenty of formatting options (two-digit vs. four-digit years, DMY vs MDY vs YMD, months shown as digits, time shown in decimal hours/decimal minutes/decimal seconds/decimal days, etc.)





NUTATION.CPP: Computes the nutation of the earth's equator.





OBLIQUIT.CPP: Computes the obliquity of the earth's equator for dates from -8000 to +12000.





PERSIAN.CPP: Very special-purpose code, written to demonstrate how the Persian calendar is computed. Most people will either skip this, or just use the functions in DATE.CPP without really needing to know how the functions were put together.





PLUTO.CPP: Code to compute the position of Pluto, using the series given in Meeus' Astronomical Algorithms. This uses the same data file (VSOP.BIN) as the VSOPSON.CPP code discussed below.





PRECESS.CPP: Code to build the matrix used to convert coordinates between epochs, plus code to use those matrices to convert vectors, RA/dec coordinates, and ecliptic coordinates between epochs.





REFRACT.CPP: Contains functions to convert an observed altitude (one affected by refraction) to a true altitude (one that would be seen on an airless planet), or vice versa. Two different methods are provided. (See next paragraph for a third method.)





REFRACT4.CPP: Contains a more complex/sophisticated method of computing refraction, as given in the Explanatory Supplement to the Astronomical Almanac . This one does a numerical integration through the atmosphere, including effects due to observer altitude, temperature, barometric pressure, humidity, and wavelength of light being refracted. (But the end result is usually not all that different from that given by the two methods in the previous paragraph, except for very low altitudes or extreme conditions.)





RELATIVI.CPP: Code to test out the precession of Mercury's orbit due to general relativity, and also to test a quick and easy way to include first-order general relativity in orbit computations. The code may also be of interest because it integrates using the method of Encke, with a pretty good Runge-Kutta numerical integrator. Unfortunately, it could stand some documenting...





ROCKS.CPP: Code to compute positions for 25 faint, inner satellites of the gas giants. These are all modelled as precessing ellipses.





SHOWELEM.CPP: Code to take a set of orbital elements, and produce a human-readable text version, mostly resembling the MPC's "standard" eight-line format.





SSATS.CPP: Code to compute the positions of eight of Saturn's satellites (Mimas, Enceladus, Tethys, Dione, Rhea, Titan, Hyperion, Japetus), using the theory due to G. Dourneau.





VISLIMIT.CPP: This code computes the limiting visual magnitude, sky brightness, and extinction coefficients. It's basically lifted/ported from Brad Schaefer's article and code on pages 57-60, May 1998 Sky %26amp; Telescope, "To the Visual Limits".





VSOPSON.CPP: Code to compute planetary positions from a truncated VSOP series. That series is stored in binary form on the Guide CD-ROM, and you have to load it into a buffer to use this function. The logic is very similar to that used for the full VSOP (see BIGVSOP.CPP), except that that function works by reading in the file of coefficients on an "as needed" basis. This function (and several others listed on this page) require that you download this file of coefficient data.








Code to access astronomical datasets





UCAC-2 access: Significant enough that there is a separate page for the UCAC-2 code.





CMC-14 access: Significant enough that there is a separate page for the CMC-14 code.





Ax.0 access: Source code described and provided here.





USNO-B1.0 access: Some not well-documented (but probably still helpful) code to access USNO-B1.0 in its "original" 80GByte form is available here.





Tycho-2 CD access: Given the delays in getting Guide 8 out the door, I've had requests to provide some means of accessing the Tycho-2 CD-ROMs. (These disks have the additional advantage of being free, provided by ESA and USNO. Click here for details. Though I wouldn't be surprised if they're running low on disks.)





While I haven't done much with this data (except process it for use in Guide 8, of course), I have written a small utility to access Tycho-2 data for a given star. Click here to download the DOS software, source code, and documentation. Given this, you can extract positions, proper motion, and magnitudes (both Tycho and Johnson), with error estimates, for a given Tycho-2 star. (Both the transformation to the Johnson system and the error estimates are important. "Raw" Tycho mags match Johnson magnitudes only in a loose manner, and the Tycho-2 magnitude data runs the full range from extremely accurate to moderately accurate... you _must_ check error estimates to determine which is which.)





Artificial satellites: This has grown to the point where it deserved a page of its own. Click here for details on C/C++ code to compute artificial satellite ephemerides. (The currently-posted version is a nice improvement over some older artificial satellite code I posted previously. I'm leaving the old stuff because it's still getting a little use here and there.)





Symplectic integration: This single source file is derived from code posted by David Whysong on the newsgroup sci.astro.amateur. I simplified the code a bit and added a "test code" main() function. You can now get the results of our collaboration at David's Web site. It's a very interesting way of doing numerical integration, on several counts.





The major use for it is in cases where one wants to maintain certain constants of motion (total energy, angular momentum, etc.) Usual methods of numerical integration (Runge-Kutta, Adams-Bashforth, etc.) don't necessarily do a good job of this over really long integration spans. This was the reason for David's interest in the method. He is a graduate student working on globular cluster simulations. If, say, the total energy of his simulated globular increased (or decreased) over billions of years, we may assume that the whole thing would mysteriously explode with stars flying everywhere (or collapse down to a black hole). More details are given in the source code.





DIST.CPP: This code resulted from a discussion on the newsgroup sci.astro.amateur about how to compute the geodetic distance between two lat/lon points on the earth's surface. It incorporates two methods: one assumes a spherical earth, the other (much more complex) uses a method from





PAGE2.CPP: This code can take an RA/dec position and figure out the corresponding pages to be used in the Millennium Star Atlas, Uranometria, and Sky Atlas 2000. There's also a little snippet that takes a (lunar) lat/lon and figures out the corresponding page to be used from Antonín Rükl's Atlas of the Moon.





Source code for Charon: I've started the process of turning the Charon astrometry software into an open-source program. This begins by posting the source code for Charon (about 84 KBytes). You will also need the basic astronomical source code discussed above.





Why the source for Charon is being posted


Current state of the Charon source code


Overall structure of Charon


Why the source for Charon is being posted: I do not now, nor have I ever, regarded astrometry software as a commercially profitable idea. By a slight margin, I'm able to keep myself fed and indoors through sales of Guide; trying to sell Charon is not apt to help very much. By posting it, I hope to accomplish the following:





-- People will be able to improve/modify the source code to accomplish tasks of their own. (The basic idea behind almost all open source, of course.) With any luck, some of that code will return to me and I'll post it here for others to use. I know that some people have written Charon-like code and functions for their own project; I'm hoping that we can stop re-inventing wheels.





-- The number of users of Charon has risen to the point where I'm not doing too well at keeping up with requests for improvement. I've become a bottleneck to progress.





-- Posting the source code will probably encourage somebody to start work on making a version of Charon that does not require a Guide CD. If it can access catalog data from the Internet, then it would require no CD-ROM at all. In this case, people could do astrometry with Totally Free Software. This might encourage some people to enter the hobby.





-- Verifying the fact that the code does what I tell people it does becomes easier. Some people build their own astrometry software because they would rather not have to place trust in me and my software. Such trust would be well-placed, especially since the results of the software have been examined by dozens if not hundreds of people. (The worst risk in 'rolling your own' astrometry software is that it will usually be tested only by you, and not by dozens if not hundreds of people.) But being able to check the source code ought to relieve some anxieties.





Current state of the Charon source code: All of the rest of the source posted on this site has been cleaned up nicely and documented well enough to make it useful to somebody. Right now, this is not in the least true of Charon. You can use the source to satisfy curiosity about parts of Charon, and maybe to grab bits and pieces from the code. But there are some things I need to do to make it usable.





Charon is compiled with WATCOM C/C++, the only compiler I've got that can produce 32-bit DOS code that can access graphics. (Micro$oft's compiler can produce 32-bit DOS console applications, but none that access graphics.) I'm in the process of getting the make file into a usable form; right now, it runs across assorted directories and is not 'ready to go out of the box'.





There are parts of the code that run for hundreds of lines without the faintest trace of the residue of a comment. The code isn't particularly ugly (except for the enormous main() and parts of matcher.cpp), but figuring out what does what is a frightening prospect right now.





Overall structure of Charon: The program is overdue to be broken up into certain components. To a minor extent, this breakup is present now, but it's incomplete and Charon is really one big monolith at heart. This is really a shame, because astrometry software lends itself naturally to being broken into a set of component tasks.





Breaking Charon into components would make following the structure of the code much easier. It would let us rip out certain parts, replace them, or use them by themselves. The components (or, if you prefer, "near-components") are:





(1) Load an image from a file


(2) Find and build a list of stars in an image


(3) Get the RA/dec of an object at a given time


(4) Get stars from a catalog within a given radius of a given RA/dec


(5) Pattern-match the list of image stars to the list of catalog stars


(6) Display the results and interact with user


(7) Assorted minor functions


(8) (FUTURE, not yet implemented) Automated object detection


Load an image from a file: A piece that loads an image, be it FITS, SBIG, or other, into a standardized in-memory format. This is currently done using the load_image() function in load_img.cpp. A side note: it's possible that, in preference to writing one function that can load a zillion image formats, there should be a separate routine that can convert a zillion image formats into FITS. If I had it to do over again, I'd do it that way.





Find and build a list of stars in an image: Takes an image in the aforementioned standardized in-memory format and return a list of stars in it, with their pixel coordinates and brightnesses. The code in findstar.cpp does this right now. (I'm not totally thrilled with the job it does on noisy images; it often finds hordes of spurious stars in noise. Ideally, the code I wrote might get pitched in favor of SExtractor. I started playing with SExtractor and gave up, though; I had no real success in puzzling out what was supposed to do what.)





Get the RA/dec of an object at a given time: Given an object ID such as "V4191 Sgr" or "P/Linear S4" or "1997 XF11", and a date/time and observer position, this figures out the RA/dec of that object. The code in find_obj.cpp does this right now. If somebody tells Charon, for example, "Here's an image of asteroid 4179", Charon can figure out a rough RA/dec for 4179 at the time the image was taken, so it knows where to look for a pattern match.





In some ways, this is just a convenience. People can and do feed Charon images and just say, "Here's an image taken at about 14h51m12s, +10 43' 56"... go process it." But the 'convenience' is a major one.





Get stars from a catalog within a given radius of a given RA/dec: Given a catalog name and an RA/dec (normally provided by the preceding piece, but not always) and a radius, this piece grabs catalog data for that area and returns it in a standard way. Right now, assorted functions in grab_gsc.cpp do that. They can extract ACT, GSC, or GSC-ACT data from the Guide disks, or USNO Ax.0 or SAx.0 data from the CDs distributed by USNO, or GSC 1.2 data from the files that come from the STScI server. All the functions feed through the grab_catalog_stars() function in charon.cpp.





For Guide users, this is pretty nice; it means they can get all the data they need (comet/asteroid orbits, plus the GSC, plus other handy catalogs) all on one disk. For an open-source project, though, it'll be nice if the 'catalog grabbing' code can 'grab' from the two GSC CD-ROMs distributed by STScI, or via Internet. The first, at least, ought to be very straightforward. The second will be slightly trickier, but opens the program to anyone connected to the Internet... this would be a huge plus. I note that VizieR, for example, provides a way to ask for a given list of data from a catalog, specified by RA/dec and radius, using a single URL... almost exactly what we'd need for a stunt like this.





Right now, people have to shell out some money to buy the Guide disk, or the two GSC disks, to get into astrometry. Neither is really pricey, as astronomy products go ($89 and $69, respectively). But if trying out astrometry software was totally free, it might draw in some fresh blood.





Pattern-match the list of image stars to the list of catalog stars: Given a list of image stars (probably from the 'get stars from image' component), and a list of catalog stars (probably from the 'get stars from catalog' component), this part does the automated pattern-matching required to say, "Use the following transformation to go from image to RA/dec space and vice versa." The code for this is currently in matcher.cpp, and is a horror to behold... it's the sort of neighborhood where programmers only go in pairs.





Display the results and interact with user: This piece takes the image, lists of image stars and catalog stars, and transformation, and display the whole mess on the device of your choice. At present, this is done partly with code in the main( ) of charon.cpp, with a big chunk of it in dispimg.cpp and smaller pieces throughout the code. When applied to multiple images, this could allow for blinking. (The current blinking scheme in Charon is considerably less logical. I'd ignore it if possible.) This piece is intimately hooked into the entire user interface problem of how people can click on the image and manipulate it on-screen.





The previous components could be written in totally 'clean', generic ANSI source. Porting them to other OSes would be relatively straightforward. Unfortunately, the 'display/interact' component will be messy, heavily OS-dependent code.





Assorted minor functions: A whole slew of small pieces to do things such as write out MPC headers and reports, figure out which image star or catalog star is nearest to a given cursor position, provide ways to adjust settings, handle the list of MPC stations...





(FUTURE, not yet implemented) Automated object detection: Once the pattern matching code has done its work for a particular image, we can get a list of the RA/dec/magnitude values in that image. (Such lists would have value in and of themselves, of course.) Then we can take two or three such lists at different instants, and look for moving objects.


No comments:

Post a Comment