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
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
|
include <fset.h>
include <ctype.h>
include <math.h>
include <pkg/skywcs.h>
define MAX_FIELDS 100 # Maximum number of fields in list
define TABSIZE 8 # Spacing of tab stops
# T_CCTRAN -- Transform a list of x and y coordinates to RA nad DEC or vice
# versa using the celestial coordinate transformation computed by the CCMAP
# task.
procedure t_cctran()
bool forward
int inlist, outlist, reclist, geometry, xcolumn, ycolumn, min_sigdigits
int infd, outfd, lngunits, latunits
pointer sp, infile, outfile, record, xformat, yformat, str, dt
pointer sx1, sy1, sx2, sy2, coo, mw
bool clgetb(), streq()
int fntopnb(), imtopenp(), fntlenb(), fntgfnb(), clgwrd(), clgeti()
int open(), imtgetim (), imtlen()
pointer dtmap()
begin
# Allocate memory for transformation parameters structure
call smark (sp)
call salloc (infile, SZ_FNAME, TY_CHAR)
call salloc (outfile, SZ_FNAME, TY_CHAR)
call salloc (record, SZ_FNAME, TY_CHAR)
call salloc (xformat, SZ_FNAME, TY_CHAR)
call salloc (yformat, SZ_FNAME, TY_CHAR)
call salloc (str, SZ_LINE, TY_CHAR)
# Open the input and output file lists.
call clgstr ("input", Memc[str], SZ_FNAME)
if (Memc[str] == EOS)
call strcpy ("STDIN", Memc[str], SZ_FNAME)
inlist = fntopnb(Memc[str], NO)
call clgstr ("output", Memc[str], SZ_FNAME)
if (Memc[str] == EOS)
call strcpy ("STDOUT", Memc[str], SZ_FNAME)
outlist = fntopnb (Memc[str], NO)
call clgstr ("database", Memc[str], SZ_FNAME)
if (Memc[str] != EOS) {
dt = dtmap (Memc[str], READ_ONLY)
reclist = imtopenp ("solution")
} else {
dt = NULL
reclist = NULL
}
# Test the input and out file and record lists for validity.
if (fntlenb(inlist) <= 0)
call error (0, "The input file list is empty")
if (fntlenb(outlist) <= 0)
call error (0, "The output file list is empty")
if (fntlenb(outlist) > 1 && fntlenb(outlist) != fntlenb(inlist))
call error (0,
"Input and output file lists are not the same length")
if (dt != NULL && reclist != NULL) {
if (imtlen (reclist) > 1 && imtlen (reclist) != fntlenb (inlist))
call error (0,
"Input file and record lists are not the same length.")
}
# Get the fitting geometry.
geometry = clgwrd ("geometry", Memc[str], SZ_LINE,
",linear,distortion,geometric,")
forward = clgetb ("forward")
# Get the input and output file parameters.
iferr (lngunits = clgwrd ("lngunits", Memc[str], SZ_LINE,
SKY_LNG_UNITLIST))
lngunits = 0
iferr (latunits = clgwrd ("latunits", Memc[str], SZ_LINE,
SKY_LAT_UNITLIST))
latunits = 0
xcolumn = clgeti ("xcolumn")
ycolumn = clgeti ("ycolumn")
call clgstr ("lngformat", Memc[xformat], SZ_FNAME)
call clgstr ("latformat", Memc[yformat], SZ_FNAME)
min_sigdigits = clgeti ("min_sigdigits")
# Get the output file name.
if (fntgfnb (outlist, Memc[outfile], SZ_FNAME) == EOF)
call strcpy ("STDOUT", Memc[outfile], SZ_FNAME)
outfd = open (Memc[outfile], NEW_FILE, TEXT_FILE)
if (streq (Memc[outfile], "STDOUT") || outfd == STDOUT)
call fseti (outfd, F_FLUSHNL, YES)
# Get the record name.
if (reclist == NULL)
Memc[record] = EOS
else if (imtgetim (reclist, Memc[record], SZ_FNAME) == EOF)
Memc[record] = EOS
# Call procedure to get parameters and fill structure.
coo = NULL; sx1 = NULL; sy1 = NULL; sx2 = NULL; sy2 = NULL
call cc_init_transform (dt, Memc[record], geometry, lngunits,
latunits, sx1, sy1, sx2, sy2, mw, coo)
# While input list is not depleted, open file and transform list.
while (fntgfnb (inlist, Memc[infile], SZ_FNAME) != EOF) {
infd = open (Memc[infile], READ_ONLY, TEXT_FILE)
# Transform the coordinates.
call cc_transform_file (infd, outfd, xcolumn, ycolumn, lngunits,
latunits, Memc[xformat], Memc[yformat], min_sigdigits, sx1,
sy1, sx2, sy2, mw, coo, forward)
# Do not get a new output file name if there is not output
# file list or if only one output file was specified.
# Otherwise fetch the new name.
if (fntlenb(outlist) > 1) {
call close (outfd)
if (fntgfnb (outlist, Memc[outfile], SZ_FNAME) != EOF)
outfd = open (Memc[outfile], NEW_FILE, TEXT_FILE)
if (streq (Memc[outfile], "STDOUT") || outfd == STDOUT)
call fseti (outfd, F_FLUSHNL, YES)
}
call close (infd)
# Do not reset the transformation if there is no record list
# or only one record is specified. Otherwise fetch the next
# record name.
if (reclist != NULL && imtlen (reclist) > 1) {
if (imtgetim (reclist, Memc[record], SZ_FNAME) != EOF) {
call cc_free_transform (sx1, sy1, sx2, sy2, mw, coo)
call cc_init_transform (dt, Memc[record], geometry,
lngunits, latunits, sx1, sy1, sx2, sy2, mw, coo)
}
}
}
# Free the surface descriptors.
call cc_free_transform (sx1, sy1, sx2, sy2, mw, coo)
# Close up file and record templates.
if (dt != NULL)
call dtunmap (dt)
call close (outfd)
call fntclsb (inlist)
call fntclsb (outlist)
if (reclist != NULL)
call imtclose (reclist)
call sfree (sp)
end
# CC_TRANSFORM_FILE -- This procedure is called once for each file
# in the input list. For each line in the input file that isn't
# blank or comment, the line is transformed. Blank and comment
# lines are output unaltered.
procedure cc_transform_file (infd, outfd, xfield, yfield, lngunits,
latunits, xformat, yformat, min_sigdigits, sx1, sy1, sx2, sy2,
mw, coo, forward)
int infd #I the input file descriptor
int outfd #I the output file descriptor
int xfield #I the x column number
int yfield #I the y column number
int lngunits #I the ra / longitude units
int latunits #I the dec / latitude units
char xformat[ARB] #I output format of the x coordinate
char yformat[ARB] #I output format of the y coordinate
int min_sigdigits #I the minimum number of digits to be output
pointer sx1, sy1 #I pointers to the linear x and y surfaces
pointer sx2, sy2 #I pointers to the x and y distortion surfaces
pointer mw #I pointer to the mwcs structure
pointer coo #I pointer to the celestial coordinate structure
bool forward #I forwards transform ?
double xd, yd, xtd, ytd
int max_fields, nline, nfields, nchars, nsdig_x, nsdig_y, offset
int tlngunits, tlatunits
pointer sp, inbuf, linebuf, field_pos, outbuf, ip, ct, txformat, tyformat
int getline(), li_get_numd(), sk_stati()
pointer mw_sctran()
begin
# Allocate some working space.
call smark (sp)
call salloc (inbuf, SZ_LINE, TY_CHAR)
call salloc (linebuf, SZ_LINE, TY_CHAR)
call salloc (field_pos, MAX_FIELDS, TY_INT)
call salloc (outbuf, SZ_LINE, TY_CHAR)
call salloc (txformat, SZ_LINE, TY_CHAR)
call salloc (tyformat, SZ_LINE, TY_CHAR)
# Determine the units.
if (lngunits <= 0) {
if (coo == NULL)
tlngunits = SKY_HOURS
else
tlngunits = sk_stati (coo, S_NLNGUNITS)
} else
tlngunits = lngunits
if (latunits <= 0) {
if (coo == NULL)
tlatunits = SKY_DEGREES
else
tlatunits = sk_stati (coo, S_NLATUNITS)
} else
tlatunits = latunits
# If the formats are undefined set suitable default formats.
if (xformat[1] == EOS) {
if (! forward)
call strcpy ("%10.3f", Memc[txformat], SZ_FNAME)
else {
switch (tlngunits) {
case SKY_HOURS:
call strcpy ("%12.2h", Memc[txformat], SZ_FNAME)
case SKY_DEGREES:
call strcpy ("%11.1h", Memc[txformat], SZ_FNAME)
case SKY_RADIANS:
call strcpy ("%13.7g", Memc[txformat], SZ_FNAME)
default:
call strcpy ("%10.3f", Memc[txformat], SZ_FNAME)
}
}
} else
call strcpy (xformat, Memc[txformat], SZ_FNAME)
if (yformat[1] == EOS) {
if (! forward)
call strcpy ("%10.3f", Memc[tyformat], SZ_FNAME)
else {
switch (tlatunits) {
case SKY_HOURS:
call strcpy ("%12.2h", Memc[tyformat], SZ_FNAME)
case SKY_DEGREES:
call strcpy ("%11.1h", Memc[tyformat], SZ_FNAME)
case SKY_RADIANS:
call strcpy ("%13.7g", Memc[tyformat], SZ_FNAME)
default:
call strcpy ("%10.3f", Memc[tyformat], SZ_FNAME)
}
}
} else
call strcpy (yformat, Memc[tyformat], SZ_FNAME)
# If the transformation can be represented by mwcs then compile the
# appropriate transform. Other wise use the surface fitting code
# to do the transformation.
if (mw != NULL) {
if (forward)
ct = mw_sctran (mw, "logical", "world", 03B)
else
ct = mw_sctran (mw, "world", "logical", 03B)
}
max_fields = MAX_FIELDS
for (nline=1; getline (infd, Memc[inbuf]) != EOF; nline = nline + 1) {
for (ip=inbuf; IS_WHITE(Memc[ip]); ip=ip+1)
;
if (Memc[ip] == '#') {
# Pass comment lines on to the output unchanged.
call putline (outfd, Memc[inbuf])
next
} else if (Memc[ip] == '\n' || Memc[ip] == EOS) {
# Blank lines too.
call putline (outfd, Memc[inbuf])
next
}
# If the transformation is undefined then pass the line on
# undisturbed.
if (mw == NULL) {
call putline (outfd, Memc[inbuf])
next
}
# Expand tabs into blanks, determine field offsets.
call strdetab (Memc[inbuf], Memc[linebuf], SZ_LINE, TABSIZE)
call li_find_fields (Memc[linebuf], Memi[field_pos], max_fields,
nfields)
if (xfield > nfields || yfield > nfields) {
call fstats (infd, F_FILENAME, Memc[outbuf], SZ_LINE)
call eprintf ("Not enough fields in file %s line %d\n")
call pargstr (Memc[outbuf])
call pargi (nline)
call putline (outfd, Memc[linebuf])
next
}
offset = Memi[field_pos+xfield-1]
nchars = li_get_numd (Memc[linebuf+offset-1], xd, nsdig_x)
if (nchars == 0) {
call fstats (infd, F_FILENAME, Memc[outbuf], SZ_LINE)
call eprintf ("Bad x value in file '%s' at line %d:\n")
call pargstr (Memc[outbuf])
call pargi (nline)
call putline (outfd, Memc[linebuf])
next
}
offset = Memi[field_pos+yfield-1]
nchars = li_get_numd (Memc[linebuf+offset-1], yd, nsdig_y)
if (nchars == 0) {
call fstats (infd, F_FILENAME, Memc[outbuf], SZ_LINE)
call eprintf ("Bad y value in file '%s' at line %d:\n")
call pargstr (Memc[outbuf])
call pargi (nline)
call putline (outfd, Memc[linebuf])
next
}
# Transform the coordinates.
if (! forward) {
switch (tlngunits) {
case SKY_RADIANS:
xd = RADTODEG(xd)
case SKY_HOURS:
xd = 15.0d0 * xd
default:
;
}
switch (tlatunits) {
case SKY_RADIANS:
yd = RADTODEG(yd)
case SKY_HOURS:
yd = 15.0d0 * yd
default:
;
}
}
if (sx2 != NULL || sy2 != NULL)
call cc_do_transform (xd, yd, xtd, ytd, ct, sx1, sy1,
sx2, sy2, forward)
else
call mw_c2trand (ct, xd, yd, xtd, ytd)
if (forward) {
switch (tlngunits) {
case SKY_RADIANS:
xtd = DEGTORAD(xtd)
case SKY_HOURS:
xtd = xtd / 15.0d0
default:
;
}
switch (tlatunits) {
case SKY_RADIANS:
ytd = DEGTORAD(ytd)
case SKY_HOURS:
ytd = ytd / 15.0d0
default:
;
}
}
# Format the output line.
call li_pack_lined (Memc[linebuf], Memc[outbuf], SZ_LINE,
Memi[field_pos], nfields, xfield, yfield, xtd, ytd,
Memc[txformat], Memc[tyformat], nsdig_x, nsdig_y, min_sigdigits)
call putline (outfd, Memc[outbuf])
}
if (ct != NULL)
call mw_ctfree (ct)
call sfree (sp)
end
|
