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
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
|
include "sensfunc.h"
# SF_STDS -- Get the standard observations for the specified apertures.
# If ignoring aperture set all apertures to 1.
# This routine knows the output of the task STANDARD.
procedure sf_stds (standards, aps, ignoreaps, stds, nstds)
char standards # Standard star data file
pointer aps # Aperture list
bool ignoreaps # Ignore apertures?
pointer stds # Pointer to standard observations (returned)
int nstds # Number of standard observations (returned)
int i, j, fd, beam, npts, nwaves, nalloc
real exptime, airmass, wstart, wend
real wavelength, flux, dwave, count
pointer sp, image, title, std
pointer waves, fluxes, dwaves, counts, sens, fit, wts, iwts, x, y
bool rng_elementi()
int open(), fscan(), nscan(), stridxs()
errchk open, malloc, realloc
begin
call smark (sp)
call salloc (image, SZ_STDIMAGE, TY_CHAR)
call salloc (title, SZ_STDTITLE, TY_CHAR)
# Open the standard observation data file.
fd = open (standards, READ_ONLY, TEXT_FILE)
# Read the standard observations and create a structure for each one.
# The beginning of a new star is found by a line whose first word
# begins with the character '['. Otherwise the line is interpreted
# as a data line. All unrecognized formats are skipped.
nwaves = 0
nstds = 0
while (fscan (fd) != EOF) {
call gargwrd (Memc[image], SZ_STDIMAGE)
if (Memc[image] == '[') {
call gargi (beam)
call gargi (npts)
call gargr (exptime)
call gargr (airmass)
call gargr (wstart)
call gargr (wend)
call gargstr (Memc[title], SZ_STDTITLE)
if (nscan() < 7)
next
if (!rng_elementi (aps, beam))
next
if (IS_INDEF (exptime) || exptime <= 0.) {
call eprintf (
"%s: Warning - exposure time missing or zero, using 1 second\n")
call pargstr (Memc[image])
exptime = 1.
}
# For the first one create the pointer to the array of
# structures. For the following stars increase the size
# of the pointer array and finish up the previous standard
# star.
if (nstds == 0) {
nstds = nstds + 1
call calloc (stds, nstds, TY_INT)
call calloc (std, LEN_STD, TY_STRUCT)
Memi[stds+nstds-1] = std
} else {
if (nwaves > 0) {
call realloc (waves, nwaves, TY_REAL)
call realloc (fluxes, nwaves, TY_REAL)
call realloc (dwaves, nwaves, TY_REAL)
call realloc (counts, nwaves, TY_REAL)
call realloc (wts, nwaves, TY_REAL)
call malloc (sens, nwaves, TY_REAL)
call malloc (fit, nwaves, TY_REAL)
call malloc (iwts, nwaves, TY_REAL)
call malloc (x, nwaves, TY_REAL)
call malloc (y, nwaves, TY_REAL)
call amovr (Memr[wts], Memr[iwts], nwaves)
STD_NWAVES(std) = nwaves
STD_WAVES(std) = waves
STD_FLUXES(std) = fluxes
STD_DWAVES(std) = dwaves
STD_COUNTS(std) = counts
STD_SENS(std) = sens
STD_FIT(std) = fit
STD_WTS(std) = wts
STD_IWTS(std) = iwts
STD_X(std) = x
STD_Y(std) = y
nstds = nstds + 1
call realloc (stds, nstds, TY_INT)
call calloc (std, LEN_STD, TY_STRUCT)
Memi[stds+nstds-1] = std
}
}
# Start a new standard star.
std = Memi[stds+nstds-1]
if (ignoreaps)
STD_BEAM(std) = 1
else
STD_BEAM(std) = beam
STD_NPTS(std) = npts
STD_EXPTIME(std) = exptime
STD_AIRMASS(std) = airmass
STD_WSTART(std) = wstart
STD_WEND(std) = wend
STD_SHIFT(std) = 0.
STD_NWAVES(std) = 0
# Decode the image and sky strings.
call strcpy (Memc[title], STD_TITLE(std), SZ_STDTITLE)
i = stridxs ("]", Memc[image])
if (Memc[image+i] == ']')
i = i + 1
Memc[image+i-1] = EOS
call strcpy (Memc[image+1], STD_IMAGE(std), SZ_STDIMAGE)
if (Memc[image+i] == '-') {
i = i + 2
j = stridxs ("]", Memc[image+i]) + i
Memc[image+j-1] = EOS
call strcpy (Memc[image+i], STD_SKY(std), SZ_STDIMAGE)
} else
STD_SKY(std) = EOS
nwaves = 0
# Interprete the line as standard star wavelength point.
} else if (nstds > 0) {
call reset_scan()
call gargr (wavelength)
call gargr (flux)
call gargr (dwave)
call gargr (count)
if (nscan() < 3)
next
if (wavelength < min (wstart, wend) ||
wavelength > max (wstart, wend) ||
flux<=0. || dwave<=0. || count<=0.)
next
if (!rng_elementi (aps, beam))
next
nwaves = nwaves + 1
# Allocate in blocks to minimize the number of reallocations.
if (nwaves == 1) {
nalloc = 100
call malloc (waves, nalloc, TY_REAL)
call malloc (fluxes, nalloc, TY_REAL)
call malloc (dwaves, nalloc, TY_REAL)
call malloc (counts, nalloc, TY_REAL)
call malloc (wts, nalloc, TY_REAL)
} else if (nwaves > nalloc) {
nalloc = nalloc + 100
call realloc (waves, nalloc, TY_REAL)
call realloc (fluxes, nalloc, TY_REAL)
call realloc (dwaves, nalloc, TY_REAL)
call realloc (counts, nalloc, TY_REAL)
call realloc (wts, nalloc, TY_REAL)
}
# Record the data and compute the sensitivity.
Memr[waves+nwaves-1] = wavelength
Memr[fluxes+nwaves-1] = flux
Memr[dwaves+nwaves-1] = dwave
Memr[counts+nwaves-1] = count
Memr[wts+nwaves-1] = 1.
}
}
# Finish up the last standard star and close the file.
if (nstds > 0) {
STD_NWAVES(std) = nwaves
if (nwaves > 0) {
call realloc (waves, nwaves, TY_REAL)
call realloc (fluxes, nwaves, TY_REAL)
call realloc (dwaves, nwaves, TY_REAL)
call realloc (counts, nwaves, TY_REAL)
call realloc (wts, nwaves, TY_REAL)
call malloc (sens, nwaves, TY_REAL)
call malloc (fit, nwaves, TY_REAL)
call malloc (iwts, nwaves, TY_REAL)
call malloc (x, nwaves, TY_REAL)
call malloc (y, nwaves, TY_REAL)
call amovr (Memr[wts], Memr[iwts], nwaves)
STD_WAVES(std) = waves
STD_FLUXES(std) = fluxes
STD_DWAVES(std) = dwaves
STD_COUNTS(std) = counts
STD_SENS(std) = sens
STD_FIT(std) = fit
STD_WTS(std) = wts
STD_IWTS(std) = iwts
STD_X(std) = x
STD_Y(std) = y
}
}
call close (fd)
call sfree (sp)
# Add standard stars for any added points and composite points.
nstds = nstds + 2
call realloc (stds, nstds, TY_INT)
call calloc (std, LEN_STD, TY_STRUCT)
Memi[stds+nstds-2] = std
call strcpy ("Added", STD_IMAGE(std), SZ_STDIMAGE)
STD_BEAM(std) = STD_BEAM(Memi[stds])
STD_NPTS(std) = STD_NPTS(Memi[stds])
STD_EXPTIME(std) = 1.
STD_AIRMASS(std) = 1.
STD_WSTART(std) = STD_WSTART(Memi[stds])
STD_WEND(std) = STD_WEND(Memi[stds])
STD_SHIFT(std) = 0.
STD_NWAVES(std) = 0
call calloc (std, LEN_STD, TY_STRUCT)
Memi[stds+nstds-1] = std
call strcpy ("Composite", STD_IMAGE(std), SZ_STDIMAGE)
STD_BEAM(std) = STD_BEAM(Memi[stds])
STD_NPTS(std) = STD_NPTS(Memi[stds])
STD_EXPTIME(std) = 1.
STD_AIRMASS(std) = 1.
STD_WSTART(std) = STD_WSTART(Memi[stds])
STD_WEND(std) = STD_WEND(Memi[stds])
STD_SHIFT(std) = 0.
STD_NWAVES(std) = 0
end
# SF_FREE -- Free the standard observations and aperture array.
procedure sf_free (stds, nstds, apertures, napertures)
pointer stds # Pointer to standard observations
int nstds # Number of standard observations
pointer apertures # Pointer to apertures array
int napertures # Number of apertures
int i
pointer std
begin
do i = 1, nstds {
std = Memi[stds+i-1]
if (STD_NWAVES(std) > 0) {
call mfree (STD_WAVES(std), TY_REAL)
call mfree (STD_FLUXES(std), TY_REAL)
call mfree (STD_DWAVES(std), TY_REAL)
call mfree (STD_COUNTS(std), TY_REAL)
call mfree (STD_SENS(std), TY_REAL)
call mfree (STD_FIT(std), TY_REAL)
call mfree (STD_WTS(std), TY_REAL)
call mfree (STD_IWTS(std), TY_REAL)
call mfree (STD_X(std), TY_REAL)
call mfree (STD_Y(std), TY_REAL)
}
call mfree (std, TY_STRUCT)
}
call mfree (stds, TY_INT)
call mfree (apertures, TY_INT)
end
|