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
|
# SKYGROUP -- Group coordinate list on the sky.
#
# The input is list of ra (0h-24h or 0d-360d) and dec (-90d to 90d) in the
# first two columns followed by arbitrary data (usually a filename).
# This is complicated by the periodicities at 0h.
procedure skygroup (input, output)
file input {prompt="Input list"}
string output {prompt="Output rootname"}
string extn = "" {prompt="Optional output extension"}
real sep = 60 {prompt="Separation between groups (arcsec)"}
string raunit = "hr" {prompt="RA unit (hr|deg)"}
bool keepcoords = yes {prompt="Keep coordinates in output?"}
string raformat = "%.2h" {prompt="Output RA format"}
string decformat = "%.1h" {prompt="Output DEC format"}
struct *fd1, *fd2
begin
file in, out, fname, temp1, temp2, temp3
int i, j, n, n1
real dra, r1, d1, r2, d2, r3, d3, r4
struct fmtstr, data1, data2, data3
# Temporary files.
fname = mktemp ("tmp")
in = fname // "in"
temp1 = fname // "1"
temp2 = fname // "2"
temp3 = fname // "3"
# Set parameters.
fname = input
out = output
# Check for existing output files.
files (out//"_[0-9][0-9][0-9][0-9]*"//extn, > temp1)
count (temp1) | scan (n); delete (temp1, verify-)
if (access(out) || n > 0)
error (1, "Output files already exist")
if (raunit == "hr")
dra = 24
else
dra = 360
if (keepcoords)
fmtstr = "%d " // raformat // " " // decformat // " %s"
else
fmtstr = "%d %s"
fmtstr += "\n"
# We start by sorting in dec.
sort (fname, col=2, num+, > in)
# Find jumps in dec bigger than the separation and then sort
# in ra and find jumps in ra bigger than separation. Handle
# the wrap around at 0h by duplicating to extend beyond 24h.
# The duplicates will be eliminated during the merging process.
n = 1
fd1 = in
if (fscan (fd1, r1, d1, data1) == EOF)
error (1, "No data or badly formated data")
while (fscan (fd1, r2, d2, data2) != EOF) {
print (r1, d1, data1, >> temp1)
if (r1 / max (0.001, dcos(d1)) * 3600 <= sep) {
r4 = r1 + dra
print (r4, d1, data1, >> temp1)
}
if (abs(d2-d1) <= sep) {
r1 = r2; d1 = d2; data1 = data2
next
}
r3 = r2; d3 = d2; data3 = data2
sort (temp1, col=1, num+, > temp2)
delete (temp1, verify-)
fd2 = temp2
if (fscan (fd2, r1, d1, data1) == EOF);
while (fscan (fd2, r2, d2, data2) != EOF) {
if (keepcoords)
printf (fmtstr, n, r1, d1, data1, >> temp3)
else
printf (fmtstr, n, data1, >> temp3)
skysep (r1, d1, r2, d2, raunit=raunit, verb-)
if (skysep.sep <= sep) {
r1 = r2; d1 = d2; data1 = data2
next
}
n += 1
r1 = r2; d1 = d2; data1 = data2
}
fd2 = ""; delete (temp2, verify-)
if (keepcoords)
printf (fmtstr, n, r1, d1, data1, >> temp3)
else
printf (fmtstr, n, data1, >> temp3)
n += 1
r1 = r3; d1 = d3; data1 = data3
}
fd1 = ""; delete (in, verify-)
print (r1, d1, data1, >> temp1)
if (r1 / max (0.001, dcos(d1)) * 3600 <= sep) {
r4 = r1 + dra
print (r4, d1, data1, >> temp1)
}
sort (temp1, col=1, num+, > temp2)
delete (temp1, verify-)
fd2 = temp2
if (fscan (fd2, r1, d1, data1) == EOF);
while (fscan (fd2, r2, d2, data2) != EOF) {
if (keepcoords)
printf (fmtstr, n, r1, d1, data1, >> temp3)
else
printf (fmtstr, n, data1, >> temp3)
skysep (r1, d1, r2, d2, raunit=raunit, verb-)
if (skysep.sep <= sep) {
r1 = r2; d1 = d2; data1 = data2
next
}
n += 1
r1 = r2; d1 = d2; data1 = data2
}
fd2 = ""; delete (temp2, verify-)
if (keepcoords)
printf (fmtstr, n, r1, d1, data1, >> temp3)
else
printf (fmtstr, n, data1, >> temp3)
# Now write out the lists and check for duplicate which must be
# merged.
sort (temp3, col=1, num+, > temp1); delete (temp3, verify-)
touch (temp2)
fd1 = temp1
if (fscan (fd1, i, data1) == EOF);
while (fscan (fd1, j, data2) != EOF) {
if (data1 == data2) {
print (j, i, >> temp2)
next
}
printf ("%s_%03d%s\n", out, i, extn) | scan (fname)
print (data1, >> fname)
i = j; data1 = data2
}
fd1 = ""; delete (temp1, verify-)
printf ("%s_%03d%s\n", out, i, extn) | scan (fname)
print (data1, >> fname)
sort (temp2, col=1, num+, rev+) | unique (> temp1)
delete (temp2, verify-)
# Merge the lists.
n1 = n
fd1 = temp1
while (fscan (fd1, j, i) != EOF) {
printf ("%s_%03d%s\n", out, j, extn) | scan (fname)
if (access (fname)) {
printf ("%s_%03d%s\n", out, i, extn) | scan (in)
concat (in, fname, append+)
delete (in, verify-)
n1 -= 1
}
}
fd1 = ""; delete (temp1, verify-)
# Renumber if needed.
if (n1 != n) {
i = 1
for (j=1; j<=n; j+=1) {
printf ("%s_%03d%s\n", out, j, extn) | scan (fname)
if (access(fname)) {
if (i != j) {
printf ("%s_%03d%s\n", out, i, extn) | scan (in)
rename (fname, in)
}
i += 1
}
}
}
# Create the final output list of lists.
files (out//"_[0-9]*", > out//extn)
end
|
