10abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#!/usr/bin/python 20abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 30abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com''' 40abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comCopyright 2013 Google Inc. 50abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 60abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comUse of this source code is governed by a BSD-style license that can be 70abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comfound in the LICENSE file. 80abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com''' 90abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 100abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comimport math 110abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comimport pprint 120abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 130abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comdef withinStdDev(n): 140abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com """Returns the percent of samples within n std deviations of the normal.""" 150abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return math.erf(n / math.sqrt(2)) 160abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 170abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comdef withinStdDevRange(a, b): 180abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com """Returns the percent of samples within the std deviation range a, b""" 190abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if b < a: 200abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return 0; 210abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 220abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if a < 0: 230abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if b < 0: 240abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return (withinStdDev(-a) - withinStdDev(-b)) / 2; 250abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com else: 260abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return (withinStdDev(-a) + withinStdDev(b)) / 2; 270abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com else: 280abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return (withinStdDev(b) - withinStdDev(a)) / 2; 290abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 300abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 310abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#We have a bunch of smudged samples which represent the average coverage of a range. 320abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#We have a 'center' which may not line up with those samples. 330abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#From the 'center' we want to make a normal approximation where '5' sample width out we're at '3' std deviations. 340abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#The first and last samples may not be fully covered. 350abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 360abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#This is the sub-sample shift for each set of FIR coefficients (the centers of the lcds in the samples) 370abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#Each subpxl takes up 1/3 of a pixel, so they are centered at x=(i/n+1/2n), or 1/6, 3/6, 5/6 of a pixel. 380abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#Each sample takes up 1/4 of a pixel, so the results fall at (x*4)%1, or 2/3, 0, 1/3 of a sample. 390abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comsamples_per_pixel = 4 400abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comsubpxls_per_pixel = 3 410abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#sample_offsets is (frac, int) in sample units. 420abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comsample_offsets = [math.modf((float(subpxl_index)/subpxls_per_pixel + 1.0/(2.0*subpxls_per_pixel))*samples_per_pixel) for subpxl_index in range(subpxls_per_pixel)] 430abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 440abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#How many samples to consider to the left and right of the subpxl center. 450abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comsample_units_width = 5 460abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 470abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#The std deviation at sample_units_width. 480abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comstd_dev_max = 3 490abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 500abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#The target sum is in some fixed point representation. 510abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com#Values larger the 1 in fixed point simulate ink spread. 520abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comtarget_sum = 0x110 530abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 540abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.comfor sample_offset, sample_align in sample_offsets: 550abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs = [] 560abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs_rounded = [] 570abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 580abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com #We start at sample_offset - sample_units_width 590abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_sample_left = sample_offset - sample_units_width 600abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_std_dev_left = -std_dev_max 610abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 620abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com done = False 630abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com while not done: 640abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_sample_right = math.floor(current_sample_left + 1) 650abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if current_sample_right > sample_offset + sample_units_width: 660abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com done = True 670abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_sample_right = sample_offset + sample_units_width 680abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_std_dev_right = current_std_dev_left + ((current_sample_right - current_sample_left) / sample_units_width) * std_dev_max 690abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 700abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coverage = withinStdDevRange(current_std_dev_left, current_std_dev_right) 710abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs.append(coverage * target_sum) 720abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs_rounded.append(int(round(coverage * target_sum))) 730abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 740abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_sample_left = current_sample_right 750abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com current_std_dev_left = current_std_dev_right 760abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 770abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # Now we have the numbers we want, but our rounding needs to add up to target_sum. 780abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com delta = 0 790abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs_rounded_sum = sum(coeffs_rounded) 800abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if coeffs_rounded_sum > target_sum: 810abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # The coeffs add up to too much. Subtract 1 from the ones which were rounded up the most. 820abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com delta = -1 830abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 840abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if coeffs_rounded_sum < target_sum: 850abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # The coeffs add up to too little. Add 1 to the ones which were rounded down the most. 860abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com delta = 1 870abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 880abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com if delta: 890abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com print "Initial sum is 0x%0.2X, adjusting." % (coeffs_rounded_sum,) 900abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeff_diff = [(coeff_rounded - coeff) * delta 910abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com for coeff, coeff_rounded in zip(coeffs, coeffs_rounded)] 920abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 930abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com class IndexTracker: 940abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com def __init__(self, index, item): 950abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com self.index = index 960abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com self.item = item 970abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com def __lt__(self, other): 980abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return self.item < other.item 990abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com def __repr__(self): 1000abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com return "arr[%d] == %s" % (self.index, repr(self.item)) 1010abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 1020abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeff_pkg = [IndexTracker(i, diff) for i, diff in enumerate(coeff_diff)] 1030abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeff_pkg.sort() 1040abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 1050abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # num_elements_to_force_round had better be < (2 * sample_units_width + 1) or 1060abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # * our math was wildy wrong 1070abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # * an awful lot of the curve is out side our sample 1080abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com # either is pretty bad, and probably means the results will not be useful. 1090abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com num_elements_to_force_round = abs(coeffs_rounded_sum - target_sum) 1100abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com for i in xrange(num_elements_to_force_round): 1110abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com print "Adding %d to index %d to force round %f." % (delta, coeff_pkg[i].index, coeffs[coeff_pkg[i].index]) 1120abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs_rounded[coeff_pkg[i].index] += delta 1130abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 1140abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com print "Prepending %d 0x00 for allignment." % (sample_align,) 1150abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com coeffs_rounded_aligned = ([0] * int(sample_align)) + coeffs_rounded 1160abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com 1170abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com print ', '.join(["0x%0.2X" % coeff_rounded for coeff_rounded in coeffs_rounded_aligned]) 1180abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com print sum(coeffs), hex(sum(coeffs_rounded)) 1190abbff9987b9452fd30cce198bea34fdb210ac41bungeman@google.com print 120