diff --git a/gstlal-inspiral/bin/gstlal_inspiral_compress_ranking_stat b/gstlal-inspiral/bin/gstlal_inspiral_compress_ranking_stat index 1ef2e054bc78d3ac8b2a5dd4a9d9dc2411db7324..830b9f5f635a8642542a6581762b0e49de87489a 100755 --- a/gstlal-inspiral/bin/gstlal_inspiral_compress_ranking_stat +++ b/gstlal-inspiral/bin/gstlal_inspiral_compress_ranking_stat @@ -54,7 +54,7 @@ def parse_command_line(): parser = OptionParser( version = "Name: %%prog\n%s" % "" # FIXME ) - parser.add_option("-t", "--threshold", type = "float", default = 0.01, help = "Only keep horizon distance values that differ by this much, fractionally, from their neighbours (default = 0.01).") + parser.add_option("-t", "--threshold", type = "float", default = 0.03, help = "Only keep horizon distance values that differ by this much, fractionally, from their neighbours (default = 0.03).") parser.add_option("-v", "--verbose", action = "store_true", help = "Be verbose.") options, filenames = parser.parse_args() @@ -115,34 +115,32 @@ for filename in filenames: abs_ln_thresh = math.log1p(options.threshold) for instrument, horizon_history in list(rankingstat.numerator.horizon_history.items()): - # GPS times - keys = horizon_history.keys() - # distances - values = horizon_history.values() - - # indexes of distances that are non-zero and differ - # fractionally from both neighbours by less than the - # selected threshold. always keep the first and last - # values - removals = [i for i in range(2, len(values) - 1) if values[i - 1] > 0. and values[i] > 0. and values[i + 1] > 0. and abs(math.log(values[i] / values[i - 1])) < abs_ln_thresh and abs(math.log(values[i] / values[i + 1])) < abs_ln_thresh] - - # compress arrays - j = k = 0 - for i in range(len(keys)): - if k < len(removals) and i == removals[k]: - k += 1 - else: - keys[j] = keys[i] - values[j] = values[i] - j += 1 - assert (i + 1) - j == len(removals), "%d %d %d" % (i, j, len(removals)) - del keys[j:] - del values[j:] + # GPS time / distance pairs + items = horizon_history.items() + + # compress array + j = 1 + for i in range(1, len(items) - 1): + values = items[j - 1][1], items[i][1], items[i + 1][1] + # remove distances that are non-zero and differ + # fractionally from both neighbours by less than + # the selected threshold. always keep the first + # and last values + if values[0] > 0. and values[1] > 0. and values[2] > 0. and abs(math.log(values[1] / values[0])) < abs_ln_thresh and abs(math.log(values[1] / values[2])) < abs_ln_thresh: + continue + # remove distances that are 0 and surrounded by 0 + # on both sides (basically the same as the last + # test, but we can't take log(0)). + if values == (0., 0., 0.): + continue + items[j] = items[i] + j += 1 + del items[j:] if options.verbose: print >>sys.stderr, "\"%s\": %s horizon history reduced to %.3g%% of original size" % (filename, instrument, 100. * j / (i + 1.)) # replace - rankingstat.numerator.horizon_history[instrument] = type(horizon_history)(zip(keys, values)) + rankingstat.numerator.horizon_history[instrument] = type(horizon_history)(items) # # re-insert into XML tree