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
|
# Copyright(c) 1986 Association of Universities for Research in Astronomy Inc.
include "../icombine.h"
# IC_MEDIAN -- Median of lines
procedure ic_medianr (d, n, npts, median)
pointer d[ARB] # Input data line pointers
int n[npts] # Number of good pixels
int npts # Number of output points per line
real median[npts] # Median
int i, j1, j2, j3, k, n1
bool even
real val1, val2, val3
include "../icombine.com"
begin
if (dflag == D_NONE) {
do i = 1, npts
median[i]= blank
return
}
# Check for previous sorting
if (mclip) {
if (dflag == D_ALL) {
n1 = n[1]
even = (mod (n1, 2) == 0)
j1 = n1 / 2 + 1
j2 = n1 / 2
do i = 1, npts {
k = i - 1
if (even) {
val1 = Memr[d[j1]+k]
val2 = Memr[d[j2]+k]
median[i] = (val1 + val2) / 2.
} else
median[i] = Memr[d[j1]+k]
}
} else {
do i = 1, npts {
k = i - 1
n1 = n[i]
if (n1 > 0) {
j1 = n1 / 2 + 1
if (mod (n1, 2) == 0) {
j2 = n1 / 2
val1 = Memr[d[j1]+k]
val2 = Memr[d[j2]+k]
median[i] = (val1 + val2) / 2.
} else
median[i] = Memr[d[j1]+k]
} else
median[i] = blank
}
}
return
}
# Repeatedly exchange the extreme values until there are three
# or fewer pixels.
do i = 1, npts {
k = i - 1
n1 = n[i]
while (n1 > 3) {
j1 = 1
j2 = 1
val1 = Memr[d[j1]+k]
val2 = val1
do j3 = 2, n1 {
val3 = Memr[d[j3]+k]
if (val3 > val1) {
j1 = j3
val1 = val3
} else if (val3 < val2) {
j2 = j3
val2 = val3
}
}
j3 = n1 - 1
if (j1 < j3 && j2 < j3) {
Memr[d[j1]+k] = val3
Memr[d[j2]+k] = Memr[d[j3]+k]
Memr[d[j3]+k] = val1
Memr[d[n1]+k] = val2
} else if (j1 < j3) {
if (j2 == j3) {
Memr[d[j1]+k] = val3
Memr[d[n1]+k] = val1
} else {
Memr[d[j1]+k] = Memr[d[j3]+k]
Memr[d[j3]+k] = val1
}
} else if (j2 < j3) {
if (j1 == j3) {
Memr[d[j2]+k] = val3
Memr[d[n1]+k] = val2
} else {
Memr[d[j2]+k] = Memr[d[j3]+k]
Memr[d[j3]+k] = val2
}
}
n1 = n1 - 2
}
if (n1 == 3) {
val1 = Memr[d[1]+k]
val2 = Memr[d[2]+k]
val3 = Memr[d[3]+k]
if (val1 < val2) {
if (val2 < val3) # abc
median[i] = val2
else if (val1 < val3) # acb
median[i] = val3
else # cab
median[i] = val1
} else {
if (val2 > val3) # cba
median[i] = val2
else if (val1 > val3) # bca
median[i] = val3
else # bac
median[i] = val1
}
} else if (n1 == 2) {
val1 = Memr[d[1]+k]
val2 = Memr[d[2]+k]
median[i] = (val1 + val2) / 2
} else if (n1 == 1)
median[i] = Memr[d[1]+k]
else
median[i] = blank
}
end
|
