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
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
|
include <mach.h>
include <imhdr.h>
include "../lib/daophotdef.h"
include "../lib/psfdef.h"
# DP_PSUBRAST -- Fetch the prospective PSF star data and check it for bad
# pixel values.
pointer procedure dp_psubrast (dao, im, lowbad, highbad, x1, x2, y1, y2,
saturated)
pointer dao # pointer to the daophot structure
pointer im # pointer to the input image
real lowbad, highbad # minimum and maximum good data values
int x1, x2 # output x limits of the extracted subraster
int y1, y2 # output y limits of the extracted subraster
int saturated # is the star saturated ?
pointer psf, buf
real psfrad, fitrad
int dp_chksr()
pointer dp_subrast()
begin
# Initialize.
psf = DP_PSF(dao)
buf = NULL
psfrad = DP_PSFRAD(dao)
fitrad = DP_FITRAD(dao)
# Get the data.
buf = dp_subrast (im, DP_CUR_PSFX(psf), DP_CUR_PSFY(psf), psfrad,
x1, x2, y1, y2)
if (buf == NULL)
return (NULL)
# Check for bad pixels in subraster, compute the min and max.
if (dp_chksr (DP_CUR_PSFX(psf), DP_CUR_PSFY(psf), Memr[buf],
x2 - x1 + 1, y2 - y1 + 1, x1, y1, psfrad, fitrad,
lowbad, highbad, saturated, DP_CUR_PSFMIN(psf),
DP_CUR_PSFMAX(psf), DP_CUR_PSFGMAX(psf)) == ERR) {
call mfree (buf, TY_REAL)
return (NULL)
}
return (buf)
end
# DL_LSUBRAST -- Give a valid PSF star compute the limits of the data to
# be extracted around it.
int procedure dp_lsubrast (im, xcen, ycen, radius, x1, x2, y1, y2)
pointer im # input image descriptor
real xcen, ycen # center of subraster
real radius # radius of the box
int x1, y1, x2, y2 # boundaries of subraster
begin
# Calculate start position of extraction box.
x1 = int (xcen - radius) - 2
x2 = int (xcen + radius) + 3
y1 = int (ycen - radius) - 2
y2 = int (ycen + radius) + 3
if (x1 > IM_LEN(im,1) || x2 < 1 || y1 > IM_LEN(im,2) || y2 < 1)
return (ERR)
x1 = max (1, x1)
x2 = min (IM_LEN(im,1), x2)
y1 = max (1, y1)
y2 = min (IM_LEN(im,2), y2)
return (OK)
end
# DP_SUBRAST -- Given a valid PSF star extract the data around it.
pointer procedure dp_subrast (im, xcen, ycen, radius, x1, x2, y1, y2)
pointer im # input image descriptor
real xcen, ycen # center of subraster
real radius # radius of the box
int x1, y1, x2, y2 # boundaries of subraster
int j, ncols
pointer buf, ptr, imbuf
pointer imgs2r()
begin
# Calculate start position of extraction box.
x1 = int (xcen - radius) - 2
x2 = int (xcen + radius) + 3
y1 = int (ycen - radius) - 2
y2 = int (ycen + radius) + 3
if (x1 > IM_LEN(im,1) || x2 < 1 || y1 > IM_LEN(im,2) || y2 < 1)
return (NULL)
x1 = max (1, x1)
x2 = min (IM_LEN(im,1), x2)
y1 = max (1, y1)
y2 = min (IM_LEN(im,2), y2)
call malloc (buf, (x2 - x1 + 1) * (y2 - y1 + 1), TY_REAL)
ptr = buf
ncols = x2 - x1 + 1
do j = y1, y2 {
imbuf = imgs2r (im, x1, x2, j, j)
call amovr (Memr[imbuf], Memr[ptr], ncols)
ptr = ptr + ncols
}
return (buf)
end
# DP_CHKSR -- Check the input subraster for bad pixels.
int procedure dp_chksr (x, y, sr, xdim, ydim, x1, y1, psfrad, fitrad, lowbad,
highbad, saturated, dmin, dmax, gmax)
real x, y # position of the star
real sr[xdim,ydim] # the data subraster
int xdim, ydim # the dimensions of the subraster
int x1, y1 # the lower left hand coordinates of the array
real psfrad # the psf radius
real fitrad # the fitting radius
real lowbad, highbad # the good data limits
int saturated # is the star saturated
real dmin, dmax # output data limits
real gmax # maximum good data limit
int i,j
real pradsq, fradsq, dy2, r2
begin
pradsq = psfrad * psfrad
fradsq = fitrad * fitrad
dmin = MAX_REAL
dmax = -MAX_REAL
gmax = -MAX_REAL
saturated = NO
# Loop through the pixels looking for bad values.
do j = 1, ydim {
dy2 = (y - (y1 + j -1)) ** 2
if (dy2 > pradsq)
next
do i = 1, xdim {
r2 = (x - (x1 + i - 1)) ** 2 + dy2
if (r2 > pradsq)
next
if (sr[i,j] < dmin)
dmin = sr[i,j]
if (sr[i,j] > dmax)
dmax = sr[i,j]
if (sr[i,j] < lowbad || sr[i,j] > highbad) {
if (r2 <= fradsq) {
if (sr[i,j] < lowbad)
return (ERR)
else if (saturated == NO)
saturated = YES
}
} else {
if (sr[i,j] > gmax)
gmax = sr[i,j]
}
}
}
return (OK)
end
|