-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathimg_scale.py
225 lines (185 loc) · 6.29 KB
/
img_scale.py
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
#
# Written by Min-Su Shin
# Department of Astrophysical Sciences, Princeton University
#
# You can freely use the code.
#
import numpy
import math
import stats
import numpy as np
def sky_median_sig_clip(input_arr, sig_fract, percent_fract, max_iter=100):
"""Estimating sky value for a given number of iterations
@type input_arr: numpy array
@param input_arr: image data array
@type sig_fract: float
@param sig_fract: fraction of sigma clipping
@type percent_fract: float
@param percent_fract: convergence fraction
@type max_iter: max. of iterations
@rtype: tuple
@return: (sky value, number of iteration)
"""
work_arr = numpy.ravel(input_arr)
old_sky = numpy.median(work_arr)
sig = work_arr.std()
upper_limit = old_sky + sig_fract * sig
lower_limit = old_sky - sig_fract * sig
indices = numpy.where((work_arr < upper_limit) & (work_arr > lower_limit))
work_arr = work_arr[indices]
new_sky = numpy.median(work_arr)
iteration = 0
while ((math.fabs(old_sky - new_sky)/new_sky) > percent_fract) and (iteration < max_iter) :
iteration += 1
old_sky = new_sky
sig = work_arr.std()
upper_limit = old_sky + sig_fract * sig
lower_limit = old_sky - sig_fract * sig
indices = numpy.where((work_arr < upper_limit) & (work_arr > lower_limit))
work_arr = work_arr[indices]
new_sky = numpy.median(work_arr)
return (new_sky, iteration)
def sky_mean_sig_clip(input_arr, sig_fract, percent_fract, max_iter=100):
"""Estimating sky value for a given number of iterations
@type input_arr: numpy array
@param input_arr: image data array
@type sig_fract: float
@param sig_fract: fraction of sigma clipping
@type percent_fract: float
@param percent_fract: convergence fraction
@type max_iter: max. of iterations
@rtype: tuple
@return: (sky value, number of iteration)
"""
work_arr = numpy.ravel(input_arr)
old_sky = numpy.mean(work_arr)
sig = work_arr.std()
upper_limit = old_sky + sig_fract * sig
lower_limit = old_sky - sig_fract * sig
indices = numpy.where((work_arr < upper_limit) & (work_arr > lower_limit))
work_arr = work_arr[indices]
new_sky = numpy.mean(work_arr)
iteration = 0
while ((math.fabs(old_sky - new_sky)/new_sky) > percent_fract) and (iteration < max_iter) :
iteration += 1
old_sky = new_sky
sig = work_arr.std()
upper_limit = old_sky + sig_fract * sig
lower_limit = old_sky - sig_fract * sig
indices = numpy.where((work_arr < upper_limit) & (work_arr > lower_limit))
work_arr = work_arr[indices]
new_sky = numpy.mean(work_arr)
return (new_sky, iteration)
def linear(inputArray, scale_min=None, scale_max=None):
"""Performs linear scaling of the input numpy array.
@type inputArray: numpy array
@param inputArray: image data array
@type scale_min: float
@param scale_min: minimum data value
@type scale_max: float
@param scale_max: maximum data value
@rtype: numpy array
@return: image data array
"""
#print "img_scale : linear"
imageData=numpy.array(inputArray, copy=True)
if scale_min == None:
scale_min = imageData.min()
if scale_max == None:
scale_max = imageData.max()
imageData = imageData.clip(min=scale_min, max=scale_max)
imageData = (imageData -scale_min) / (scale_max - scale_min)
indices = numpy.where(imageData < 0)
imageData[indices] = 0.0
indices = numpy.where(imageData > 1)
imageData[indices] = 1.0
return imageData
def sqrt(inputArray, scale_min=None, scale_max=None):
"""Performs sqrt scaling of the input numpy array.
@type inputArray: numpy array
@param inputArray: image data array
@type scale_min: float
@param scale_min: minimum data value
@type scale_max: float
@param scale_max: maximum data value
@rtype: numpy array
@return: image data array
"""
#print "img_scale : sqrt"
imageData=numpy.array(inputArray, copy=True)
if scale_min == None:
scale_min = imageData.min()
if scale_max == None:
scale_max = imageData.max()
imageData = imageData.clip(min=scale_min, max=scale_max)
imageData = imageData - scale_min
indices = numpy.where(imageData < 0)
imageData[indices] = 0.0
imageData = numpy.sqrt(imageData)
imageData = imageData / math.sqrt(scale_max - scale_min)
return imageData
def log(inputArray, scale_min=None, scale_max=None):
"""Performs log10 scaling of the input numpy array.
@type inputArray: numpy array
@param inputArray: image data array
@type scale_min: float
@param scale_min: minimum data value
@type scale_max: float
@param scale_max: maximum data value
@rtype: numpy array
@return: image data array
"""
#print "img_scale : log"
imageData=numpy.array(inputArray, copy=True)
if scale_min == None:
scale_min = imageData.min()
if scale_max == None:
scale_max = imageData.max()
factor = math.log10(scale_max - scale_min)
indices0 = numpy.where(imageData < scale_min)
indices1 = numpy.where((imageData >= scale_min) & (imageData <= scale_max))
indices2 = numpy.where(imageData > scale_max)
imageData[indices0] = 0.0
imageData[indices2] = 1.0
try :
imageData[indices1] = numpy.log10(imageData[indices1])/factor
except :
print "Error on math.log10 for ", (imageData[i][j] - scale_min)
return imageData
def MAD(input_list, max, min):
data = []
for input in input_list:
data.append(input.flatten())
med = np.median(data)
sigma_MAD = 1.4826 * stats.MAD(data, med)
vmax = med + max*sigma_MAD
vmin = med - min*sigma_MAD
return vmax, vmin
def asinh(inputArray, scale_min=None, scale_max=None, non_linear=2.0):
"""Performs asinh scaling of the input numpy array.
@type inputArray: numpy array
@param inputArray: image data array
@type scale_min: float
@param scale_min: minimum data value
@type scale_max: float
@param scale_max: maximum data value
@type non_linear: float
@param non_linear: non-linearity factor
@rtype: numpy array
@return: image data array
"""
#print "img_scale : asinh"
imageData=numpy.array(inputArray, copy=True)
if scale_min == None:
scale_min = imageData.min()
if scale_max == None:
scale_max = imageData.max()
factor = numpy.arcsinh((scale_max - scale_min)/non_linear)
indices0 = numpy.where(imageData < scale_min)
indices1 = numpy.where((imageData >= scale_min) & (imageData <= scale_max))
indices2 = numpy.where(imageData > scale_max)
imageData[indices0] = 0.0
imageData[indices2] = 1.0
imageData[indices1] = numpy.arcsinh((imageData[indices1] - \
scale_min)/non_linear)/factor
return imageData