1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
|
unit unit_dss;
{calculates to RA/DEC from a DSS image pixel position}
{By han_kleijn@hnsky.org. (c) 2001, 2002, 2003
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Lesser Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with this program (see SLA_CONDITIONS); if not, write to the
* Free Software Foundation, Inc., 59 Temple Place, Suite 330,
* Boston, MA 02111-1307 USA }
interface
uses math;
var {input}
x_coeff: array[0..19] of double;{amdx1 ..20}
y_coeff: array[0..19] of double;{amdy1 ..20}
ppo_coeff: array[0..5] of double;{ppo1 ..6}
x_pixel_offset: integer; {cnpix1}
y_pixel_offset: integer; {cnpix2}
x_pixel_size : double; {xpixelSZ}
y_pixel_size : double; {ypixelSZ}
plate_ra : double; {PLTRA ..H ..M ..S}
plate_dec : double; {PLTDEC ..SN ..D ..M ..S}
dec_sign : double;
Procedure DSSPOS (xpix ,ypix : double; var xpos, ypos : double);
implementation
Procedure DSSPOS (xpix ,ypix : double; var xpos, ypos: double);//decode full astrometrical solution Digital Sky Survey image
{ Use the Digital Sky Survey polynomial solution as documented in the
README.TXT of the Digital Sky Survey by the Space Telescope Science Institute
and
DSSPOS.C of WCStools, http://tdc-www.harvard.edu/software/wcstools/ released with a GNU LESSER GENERAL PUBLIC LICENSE Version 2.1, February 1999
DSSPOS.C is based on the astrmcal.c portion of GETIMAGE by J. Doggett and the documentation (README.TXT) distributed with the Digital Sky Survey by the Space Telescope Science Institute.
Routine to determine accurate position for pixel coordinates
returns 0 if successful otherwise 1 = angle too large for projection }
{ (* Input: *)
double xpix ; (* x pixel number (RA or long without rotation) *)
double ypix ; (* y pixel number (dec or lat without rotation) *)
(* Output: *)
double *xpos ; (* Right ascension or longitude in radians *)
double *ypos ; (* Declination or latitude in radians *)}
var
x ,y ,xmm ,ymm ,xmm2 ,ymm2 ,xmm3 ,ymm3 ,x2y2,
xi, xir, eta, etar ,raoff ,ra ,dec,
ctan ,ccos : double;
const
cons2r : double= 3600*180/pi; {206264.8062470964}
twopi : double= 2*pi;
begin
//* Convert from image pixels to plate pixels */
x := xpix + x_pixel_offset -1.0+0.5; {2013 reintroduced original -1.0+0.5 factors}
y := ypix + y_pixel_offset -1.0+0.5;
{Convert from pixels to millimeters }
xmm := ( ppo_coeff[2] -x * x_pixel_size )/1000.0;
ymm := (y * y_pixel_size - ppo_coeff[5] )/1000.0;
xmm2 := xmm * xmm ;
ymm2 := ymm * ymm ;
xmm3 := xmm * xmm2 ;
ymm3 := ymm * ymm2 ;
x2y2 := xmm2 + ymm2 ;
{Compute coordinates from x,y and plate model }
xi := x_coeff[ 0] * xmm + x_coeff[ 1]* ymm +
x_coeff[ 2] + x_coeff[ 3] * xmm2 +
x_coeff[ 4] * xmm * ymm + x_coeff[ 5] * ymm2 +
x_coeff[ 6] * (x2y2 ) + x_coeff[ 7] * xmm3 +
x_coeff[ 8] * xmm2 * ymm + x_coeff[ 9] * xmm * ymm2 +
x_coeff[10] * ymm3 + x_coeff[11] * xmm *(x2y2 )+
x_coeff[12] * xmm * x2y2 * x2y2 ;
{ Ignore magnitude and color terms+ wcs->x_coeff[13]*mag + wcs->x_coeff[14]*mag*mag + wcs->x_coeff[15]*mag*mag*mag + wcs->x_coeff[16]*mag*xmm + wcs->x_coeff[17]*mag*x2y2 + wcs->x_coeff[18]*mag*xmm*x2y2 + wcs->x_coeff[19]*color; }
eta := y_coeff[ 0] * ymm + y_coeff[ 1] * xmm +
y_coeff[ 2] + y_coeff[ 3] * ymm2 +
y_coeff[ 4] * xmm *ymm + y_coeff[ 5] * xmm2 +
y_coeff[ 6] * (x2y2 ) + y_coeff[ 7] * ymm3 +
y_coeff[ 8] * ymm2 *xmm + y_coeff[ 9] * ymm * xmm2 +
y_coeff[10] * xmm3 + y_coeff[11] * ymm *(x2y2 )+
y_coeff[12] * ymm * x2y2 * x2y2 ;
{Ignore magnitude and color terms+ wcs->y_coeff[13]*mag + wcs->y_coeff[14]*mag*mag +wcs->y_coeff[15]*mag*mag*mag + wcs->y_coeff[16]*mag*ymm + wcs->y_coeff[17]*mag*x2y2) + wcs->y_coeff[18]*mag*ymm*x2y2 +wcs->y_coeff[19]*color;}
{Convert to radians }
xir := xi / cons2r ;
etar := eta / cons2r ;
(* Convert to RA and Dec *)
ctan := sin ( plate_dec )/cos(plate_dec );{tan is sin/cos}
ccos := cos ( plate_dec );
raoff := arctan2(xir / ccos ,1.0-etar *ctan );
ra := raoff + plate_ra ;
if (ra <0.0) then ra := ra +twopi ;
xpos := ra ;
dec := arctan (cos (raoff )*((etar +ctan )/(1.0-(etar *ctan ))));
ypos := dec ;
END;
end.
|
