eic-opticks/ana/decompression.py

0001 #!/usr/bin/env python
0002 #
0003 # Copyright (c) 2019 Opticks Team. All Rights Reserved.
0004 #
0005 # This file is part of Opticks
0006 # (see https://bitbucket.org/simoncblyth/opticks).
0007 #
0008 # Licensed under the Apache License, Version 2.0 (the "License");
0009 # you may not use this file except in compliance with the License.
0010 # You may obtain a copy of the License at
0011 #
0012 #   http://www.apache.org/licenses/LICENSE-2.0
0013 #
0014 # Unless required by applicable law or agreed to in writing, software
0015 # distributed under the License is distributed on an "AS IS" BASIS,
0016 # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
0017 # See the License for the specific language governing permissions and
0018 # limitations under the License.
0019 #
0020
0021 """
0022
0023 One valued inputs that land on compression bins::
0024
0025     ipdb> p aval
0026     A()sliced
0027     A([ 15.4974,  15.4974,  15.4974, ...,  15.4974,  15.4974,  15.4974])
0028     ipdb> p bval
0029     A()sliced
0030     A([ 15.4934,  15.4934,  15.4934, ...,  15.4934,  15.4934,  15.4934])
0031
0032     ipdb> p cbins[3840:3840+10]
0033     array([ 15.4692,  15.4733,  15.4773,  15.4813,  15.4853,  15.4894,  15.4934,  15.4974,  15.5014,  15.5055])
0034
0035 """
0036 import numpy as np
0037 import os, logging
0038 log = logging.getLogger(__name__)
0039
0040
0041 class Deco(object):
0042     def __init__(self, cbins, debug=True, label="", binscale=100):
0043         self.cbins = cbins
0044         self.debug = debug
0045         self.label = label
0046         self.binscale = binscale
0047
0048     @classmethod
0049     def pbins(cls, extent=5005):
0050         return np.linspace(-extent, extent, (1 << 16) - 1 )
0051
0052     def vrange(self, vals):
0053         vmin = float(min(map(lambda _:_.min(), vals)))
0054         vmax = float(max(map(lambda _:_.max(), vals)))
0055         vone = map(lambda v:np.all(v[v == v[0]]), vals)
0056         return vmin, vmax
0057
0058     def irange(self, vals):
0059
0060         """
0061
0062
0063                       |               |>>>>>>>>>>>>>>>>>>>
0064             ..........|...............|....................
0065                       |               |
0066                      vmin            vmax
0067                       |               |
0068             <<<<<<<<<<|               |
0069                      ^                 ^
0070                     imin              imax
0071
0072             In [33]: np.where(cbins < vmin)
0073             Out[33]: (array([    0,     1,     2, ..., 65467, 65468, 65469]),)
0074
0075             In [32]: np.where(cbins > vmax)
0076             Out[32]:
0077             (array([65470, 65471, 65472, 65473,
0078
0079
0080         """
0081         vmin, vmax = self.vrange(vals)
0082         assert vmax >= vmin
0083         a = vals[0]
0084         if len(vals) > 1:
0085             b = vals[1]
0086         else:
0087             b = None
0088
0089         cmin = self.cbins[0]
0090         cmax = self.cbins[-1]
0091         assert cmax >= cmin
0092
0093         if vmin < cmin:
0094             log.warning("[%s] vmin < cmin : %s < %s " % (self.label, vmin, cmin ))
0095             log.warning("a[a < cmin] %r " % a[a < cmin] )
0096             log.warning("b[b < cmin] %r " % b[b < cmin] )
0097             imin = 0
0098         else:
0099             imin = np.where(self.cbins<=vmin)[0][-1]
0100         pass
0101         if vmax > cmax:
0102             log.warning("[%s] vmax > cmax : %s > %s " % (self.label, vmax, cmax ))
0103             log.warning("a[a > cmax] %r " % a[a > cmax] )
0104             log.warning("b[b > cmax] %r " % b[b > cmax] )
0105             imax = len(self.cbins) - 1
0106         else:
0107             imax = np.where(self.cbins>=vmax)[0][0]
0108         pass
0109
0110         assert imax >= imin
0111         return imin, imax
0112
0113
0114     def bins(self, vals, mibin=10):
0115         """
0116         :param vals: list of value arrays
0117         :param mibin: minimum number of bins to grow to
0118         :return bins: array of bin edges
0119         """
0120         imin, imax = self.irange(vals)
0121         nbin = len(self.cbins)
0122         umin = imin
0123         umax = imax
0124         inc = self.binscale/10  # when inc too big cannot grow when near the edge
0125         while True:
0126             umin = max(0, umin - inc )
0127             umax = min( nbin - 1, umax + inc )
0128             bins = self.cbins[umin:umax+1][::self.binscale]
0129             lenb = len(bins)
0130             if self.debug:
0131                 log.info("umin %d umax %d lenb %d bins %r " % (umin, umax, lenb, bins))
0132             if lenb > mibin or (umin == 0 and umax == nbin - 1) :
0133                 break
0134             pass
0135         pass
0136         return bins
0137
0138
0139
0140 def decompression_bins(cbins, vals, debug=False, label="", binscale=100):
0141     """
0142     :param cbins: full range decompressed bins
0143     :param vals:
0144
0145     Compression can be considered to be a very early (primordial) binning.
0146     To avoid artifacts all subsequent binning needs to be
0147     use bins that correspond to these originals.
0148
0149     This provides a subset of full range decompression bins,
0150     corresponding to a value range.
0151
0152     ::
0153
0154         In [3]: cbins = cf.a.pbins()
0155
0156         In [4]: cbins
0157         Out[4]: array([ -24230.625 ,  -24242.3772,  -24254.1294, ..., -794375.8706, -794387.6228, -794399.375 ])
0158
0159         In [15]: np.where( cbins <=  -60254.1294 )
0160         Out[15]: (array([ 3066,  3067,  3068, ..., 65532, 65533, 65534]),)
0161
0162
0163     """
0164
0165     dc = Deco(cbins, debug=debug, label=label, binscale=binscale)
0166     return dc.bins(vals)
0167
0168
0169
0170 def test_decompression_bins_(cbins, avals, bvals, binscale=1):
0171     dbins = decompression_bins(cbins, [avals, bvals], debug=True, binscale=binscale)
0172     # assert len(dbins) > 1
0173     log.info("avals : %s " % repr(avals))
0174     log.info("bvals : %s " % repr(bvals))
0175     log.info("dbins : %s " % repr(dbins))
0176
0177 def test_decompression_bins_0():
0178     cbins = np.linspace(-300,300,10)
0179     avals = np.repeat(300,1000)
0180     bvals = np.repeat(300,1000)
0181     test_decompression_bins_(cbins, avals, bvals, binscale=1)
0182
0183 def test_decompression_bins_1():
0184     cbins = np.arange(0, 132.001, 0.004, dtype=np.float32)
0185     avals = np.repeat( cbins[3840], 100)
0186     bvals = np.repeat( cbins[3841], 100)
0187     test_decompression_bins_(cbins, avals, bvals, binscale=100)
0188
0189 def test_decompression_bins_2():
0190     cbins = np.arange(0, 132.001, 0.004, dtype=np.float32)
0191     avals = np.repeat( cbins[-1], 100)
0192     bvals = np.repeat( cbins[-1], 100)
0193     test_decompression_bins_(cbins, avals, bvals, binscale=100)
0194
0195
0196
0197
0198
0199 if __name__ == '__main__':
0200
0201     logging.basicConfig(level=logging.INFO)
0202
0203     #test_decompression_bins_0()
0204     #test_decompression_bins_1()
0205
0206     cbins = Deco.pbins(extent=5005)
0207     avals = np.repeat( cbins[-1], 100)
0208     bvals = np.repeat( cbins[-1], 100)
0209
0210     dc = Deco(cbins, debug=True)
0211
0212     bins = dc.bins([avals, bvals])
0213
0214     print("cbins\n%s" % cbins)
0215     print(avals)
0216     print(bvals)
0217     print(bins)
0218
0219
0220