;+ ; NAME: ; oc_getstar_ucac2_test ; ; PURPOSE: (one line) ; Test the UCAC2 reader routines and make several plots of found stars ; ; DESCRIPTION: ; This routine performs several simple tests using the UCAC2 reader ; routine oc_search_pos_single_ucac2(). ; ; 1) A bullseye-style plot is made showing the number of stars within ; several search radii. ; ; 2) A UCAC2 star is looked up and given in both 1950 and 2000 ; equinox *and* epoch, and the positions printed vs. Vizier's results. ; This allows for test of the proper motions in UCAC2 catalog. ; ; Note that results agree between this code and Vizier for J2000, but ; differ slightly (<~0.2 arcsec) for B1950. This is likely because ; BPRECESS and Vizier use different precession methods; this is confirmed by ; examining precession of arbitrary 2MASS stars. Therefore, BPRECESS may ; limit the accuracy of some calculations. ; ; CATEGORY: ; Occultation ; ; CALLING SEQUENCE: ; oc_getstar_ucac2_test ; ; INPUTS: ; None ; ; OPTIONAL INPUT PARAMETERS: ; None ; ; KEYWORD INPUT PARAMETERS: ; None ; ; KEYWORD OUTPUT PARAMETERS: ; None ; ; OUTPUTS: ; Makes plots to screen and prints diagnostic messages. ; ; COMMON BLOCKS: ; None ; ; SIDE EFFECTS: ; None ; ; RESTRICTIONS: ; None ; ; EXAMPLE: ; oc_getstar_ucac2_test ; ; MODIFICATION HISTORY: ; Written 1-Feb-2006 by Henry Throop, SwRI ; Modified 23-Feb-2006 by HBT. Improved documentation and formatting. ; Modified 15-Mar-2005 by HBT. Added P.M. test, vs. Vizier. ;- pro oc_getstar_ucac2_test colors x = 0.9d y = 0.9d r = 0.02 sym = 3 a=oc_search_pos_single_ucac2(x, y, r, /verbose) print, 'Radius = ' + st(r) + ', at (x,y) = (' + st(x) + ',' + st(y) + ')' print, ' HBT benchmark 1-Feb-2006: 5733 stars found, ' + $ '0.02 rad radius centred at 0.9 rad, 0.9 rad.' print, ' HBT benchmark 1-Feb-2005: 5733, 3399, 1493, ' + $ '297 at smaller radii.' print, ' HBT benchmark 1-Feb-2005: Upper edge of plot is ' + $ 'truncated due to limits of UCAC2 survey.' a=oc_search_pos_single_ucac2(x, y, r*1, /verbose) plot, a.ra, a.dec, psym=sym, $ xrange=[x-2*r, x+2*r], yrange=[y-2*r, y+2*r], $ xtitle = 'X [radians]', ytitle = 'Y [radians]' oplot, a.ra, a.dec, psym=sym, color=1 a=oc_search_pos_single_ucac2(x, y, r*0.75, /verbose) oplot, a.ra, a.dec, psym=sym, color=2 a=oc_search_pos_single_ucac2(x, y, r*0.5, /verbose) oplot, a.ra, a.dec, psym=sym, color=3 a=oc_search_pos_single_ucac2(x, y, r*0.25, /verbose) oplot, a.ra, a.dec, psym=sym, color=5 ;;; Now, test the proper motions and ET conversions ; Check J2000 coordinates. Epoch 2000, Equinox 2000. id = ['1'] id = ['35213093'] yr = cspice_jyear() d2r = !dpi / 180.d ; multiply by to convert degrees to radians r2d = 180.d / !dpi ; multiply by to convert radians to degrees d2as = 60d*60d ; multiply by to convert degrees to accsec r2as = r2d * d2as ; multiply by to convert radians to arcsec as2r = 1/r2as ; multiply by to convert arcsec to radians star = oc_getstar_ucac2(id) print print, 'Performing motion/precession tests on UCAC 35213093' print, 'Position from star.ra' print, 'RA J2000: ' + radtohms_str(star.ra) print, 'DE J2000: ' + radtodms_str(star.dec) print print, 'Now applying PM to move from catalog ET to ET=2000 yr' star = oc_star_apply_pm(star, (2000-2000)*yr) print, 'RA J2000: ' + radtohms_str(star.ra) print, 'DE J2000: ' + radtodms_str(star.dec) print print, 'Vizier: UCAC 35212093: RA J2000 = 09 58 00.597' print, 'Vizier: UCAC 35212093: DE J2000 =+09 41 40.62' print ; Check the proper motion from 2000 to 1950, using oc_star_apply, with no ; equinox change print, 'Now applying proper motion to 1950, with no precession' star = oc_getstar_ucac2(id) star = oc_star_apply_pm(star, (1950-2000)*yr) print, 'RA 2000: ' + radtohms_str(star.ra) + ' +- ' + $ radtohms_str(star.raerr) print, 'DE 2000: ' + radtodms_str(star.dec) + ' +- ' + $ radtodms_str(star.decerr) print, 'RA 2000: ' + radtohms_str(star.ra) + $ ' Epoch 1950, Equinox 2000, PM manually calculated' print, 'DE 2000: ' + radtodms_str(star.dec) + $ ' Epoch 1950, Equinox 2000, PM manually calculated' print ; Check the precession alone, with no proper motion print, 'Now applying precession to B1950, with no proper motion' star = oc_getstar_ucac2(id) bprecess, star.ra*r2d, star.dec*r2d, ra_1950_deg, dec_1950_deg print, 'RA 1950: ' + radtohms_str(ra_1950_deg*d2r) + $ ' +- ' + radtohms_str(star.raerr) print, 'DE 1950: ' + radtodms_str(dec_1950_deg*d2r) + $ ' +- ' + radtodms_str(star.decerr) ; Check the proper motion from 2000 to 1950, using manual application, with no ; equinox change star = oc_getstar_ucac2(id) star.ra = star.ra + (1950-2000)*yr*star.radot/cos(star.dec) star.dec= star.dec+ (1950-2000)*yr*star.decdot print ; Check B1950 coordinates. Epoch 1950, Equinox 1950. Can be compared to ; Vizier. Q: does the order of application matter? I assume first PM, then ; precess? It does substantially matter which way we do it in. PM must be ; defined in some coord frame -- if it's J2000, need to use J2000. *** This ; should be a new field -- the equinox of PM! print, 'Now applying proper motion and precession to 1950' star = oc_getstar_ucac2(id) star = oc_star_apply_pm(star, (1950-2000)*yr) bprecess, star.ra*r2d, star.dec*r2d, ra_1950_deg, dec_1950_deg print, 'Proper motion per year: ' + radtohms_str(star.radot*yr) print, 'Proper motion per year: ' + radtodms_str(star.decdot*yr) print, 'RA 1950: ' + radtohms_str(ra_1950_deg*d2r) + ' +- ' + $ radtohms_str(star.raerr) print, 'DE 1950: ' + radtodms_str(dec_1950_deg*d2r) + ' +- ' + $ radtodms_str(star.decerr) print print, 'Now printing Vizier catalog positions, Epoch 1950, Equinox 1950' print, 'Vizier: UCAC 35212093: RA 1950 = 09 55 20.895' print, 'Vizier: UCAC 35212093: DE 1950 =+09 56 02.41' ; Print position errors after this ; Print PM rates print, 'pmRA: ' + st(star.radot*yr*r2as*1000) + ' mas/yr' print, 'pmDE: ' + st(star.decdot*yr*r2as*1000) + ' mas/yr' stop end ;;;;;;;;;; pro other ; ** These are programming diagnostics, and not intended to be run ** ; ; We are trying to verify using Vizier that oc_getstar_ucac2(), ; Vizier, and scat are all returning the same stars. ; Below are several routines to do this. ; ; For ease of use, we start with the stars such as 1.1, 2.1, 3.1, 4.1. ; Look up 'initial guesses' of these coordinates with oc_getstar code. id_scat = ['1.1', '2.1', '3.1', '4.1'] a=oc_getstar_ucac2(id_scat) ; Look up coords w/ scat. No conversion yet -- just see that my ; 3.1 corresponds to scat's 3.1. spawn,'scat -c ucac2 -r 200 -d -t -j -n 100 ' +st(a[0].ra*r2d) +' '+ st(a[0].dec*r2d),foo &p,(foo(where(grep(foo,'0\.00$'))))[1] ; 001.000001 0.3462278 -89.6576672 14.43 13.87 13.75 15.91 18.5 -3.5 0.00 spawn,'scat -c ucac2 -r 200 -d -t -j -n 100 ' +st(a[1].ra*r2d) +' '+ st(a[1].dec*r2d),foo &p,(foo(where(grep(foo,'0\.00$'))))[1] ; 002.000001 0.1839117 -89.2438636 13.12 12.75 12.68 14.44 17.8 -9.4 0.00 spawn,'scat -c ucac2 -r 200 -d -t -j -n 100 ' +st(a[2].ra*r2d) +' '+ st(a[2].dec*r2d),foo &p,(foo(where(grep(foo,'0\.00$'))))[1] ; 003.000001 0.0079356 -88.9432128 13.60 13.03 12.92 15.51 28.6 -5.5 0.00 ; ** Looks ok. ; Convert ID's. By manually examining the file, I *believe* that 1.1->1, 2.1->876, 3.1->34546, 4.1->8117 ; Does this line successfully reproduce this? id_ucac2=oc_convertid_ucac2(SCAT=id_scat) ; Check if conversion<->reconversion keeps the same results a2=oc_getstar_ucac2(oc_convertid_ucac2(scat=id_scat)) print, a2.mag-a.mag ; 0.0000000 0.0000000 0.0000000 0.0000000 ; 0.0000000 0.0000000 0.0000000 0.0000000 ; 0.0000000 0.0000000 0.0000000 0.0000000 ; 0.0000000 0.0000000 0.0000000 0.0000000 ; ; ** Looks ok. ; Now look up coordinates. Print out the RA/DEC of these stars, and their ID's. Use Vizier to look up their ID's, and see ; if the RA/DEC match print, a.id[0]*1d, a.ra*r2d, a.dec*r2d ; 1.0000000 876.00000 3546.0000 8117.0000 ID ; 15.910000 14.440000 15.510000 15.900000 MAG[0] ; 0.34622778 0.18391167 0.0079355556 0.052879167 RA ; -89.657667 -89.243864 -88.943213 -88.132344 DEC ; ; Vizier: ; 000.3462277 -89.6576673 15.91 ; 000.1839118 -89.2438637 14.44 ; 000.0079356 -88.9432128 15.51 ; 000.0528792 -88.1323445 15.90 ; ; ** Looks ok. end