In [1]:
%matplotlib inline
from __future__ import division, absolute_import, print_function

from IPython.html.widgets import *
import numpy as np
from numpy import log
import matplotlib.pyplot as plt
from wtmm import skeletor, perform_cwt

In [2]:
def forgery(Iterations=0, Multifractal=1):
    if Multifractal:
        turns=((0.25,0.5),(0.75, 0.25))
    else:
        turns=((0.4,0.6),(0.6, 0.4))
    first_turn, second_turn = turns
    ys = [0,1]
    ts = [0,1]
    for i in range(0, Iterations + 1):
        
        j=0
        while ts[j] < 1:
            dt = ts[j+1] - ts[j] 
            dy = ys[j+1] - ys[j]
            
            ts.insert(j+1, ts[j] + first_turn[0]*dt)
            ts.insert(j+2, ts[j] + second_turn[0]*dt)
            ys.insert(j+1, ys[j] + first_turn[1]*dy)
            ys.insert(j+2, ys[j] + second_turn[1]*dy)
            
            j += 3
    
    fig, ax = plt.subplots(figsize=(8, 8), subplot_kw={'axisbg':'#EEEEEE', 'axisbelow':True})
    ax.grid(color='w', linewidth=2, linestyle='solid')
    ax.plot(ts, ys, color='b', alpha=0.4)
    ax.set_xlim(0, 1)
    ax.set_ylim(0, 1)
    
    return fig

def forgery_prices(Iterations=10, Multifractal=1):
    if Multifractal:
        turns=((0.25,0.5),(0.75, 0.25))
    else:
        turns=((0.4,0.6),(0.6, 0.4))
    first_turn, second_turn = turns
    ys = [0,1]
    ts = [0,1]
    for i in range(0, Iterations + 1):
        
        j=0
        while ts[j] < 1:
            dt = ts[j+1] - ts[j] 
            dy = ys[j+1] - ys[j]
            
            ts.insert(j+1, ts[j] + first_turn[0]*dt)
            ts.insert(j+2, ts[j] + second_turn[0]*dt)
            ys.insert(j+1, ys[j] + first_turn[1]*dy)
            ys.insert(j+2, ys[j] + second_turn[1]*dy)
            
            j += 3
    
    return np.array(ts), np.array(ys)

cartoon_t, cartoon_s = forgery_prices(5, True)
print("cartoon's signal is {} in length".format(len(cartoon_s)))
fig, ax = plt.subplots(figsize=(8, 8), subplot_kw={'axisbg':'#EEEEEE', 'axisbelow':True})
ax.grid(color='w', linewidth=2, linestyle='solid')
ax.plot(cartoon_t, cartoon_s, color='b', alpha=0.4)
ax.set_xlim(0, 1)
ax.set_ylim(0, 1)
plt.show()


cartoon's signal is 730 in length

In [ ]:
%%time
sig = cartoon_s
sig_t = cartoon_t

w_coefs = perform_cwt(sig, width_step=0.25, max_scale=100)

In [5]:
data = np.copy(w_coefs)
bifucations = skeletor(data, smallest_scale=1, proximity=9, plot=True)



In [5]:
from collections import OrderedDict

def get_distance(rc0, rc1):
    ''' 
    Get the distance between two row-column points
    '''
    r0, c0 = rc0
    r1, c1 = rc1
    return ( (c1 - c0)**2 + (r1 - r0)**2 )**0.5

    
def get_best_match(matched_lines, cur_pts, bifucs, neighborhood_threshold=1):
    
    # lines in the hood
    matches = []
    
    threshold = (cur_pts[0][0]**0.5)*(neighborhood_threshold)//1 
    if threshold < 3:
        threshold = 3
    
    # after handling a full traversal of the row domain, we can be sure that the corona are self contianed within the graph
    for (n,v), pts in bifucs.iteritems():
        
        # Again, guard against a double match with the added benefit of speeding up the loop. 
        # Things that have already been matched are also impossible to match against but would usually be processed first
        # so skip them.
        if n in matched_lines:
            continue
            
        # grab the top most point as our anchor -- what we match against.
        anchor = pts[0]
        
        # Now we walk pts and check each against the distance to the anchor. We want to track what the current minimum distance
        # because once it begins to increase, we can stop the algorithm. Further, we can default the best_distance value
        # to the distance between the anchor and the start of pts.
        
        best_distance = get_distance(cur_pts[0], anchor)
        pts_idx = None
        
        for idx, pt in enumerate(cur_pts):
            dist = get_distance(pt, anchor)
            if dist < best_distance:
                best_distance = dist
                pts_idx = idx
            elif dist > best_distance:
                break
        
        # We only want to track distances that better the last iteration's distance
        if best_distance <= threshold:
            matches.append(((n, v), pts_idx, best_distance))
            break
    
    if matches:
        # Just defaulting here to avoid a branching if-else statement
        winning_match = matches[0]
        for match in matches:
            # test the distance of each element of the list... easier than making a custom sort
            if match[2] < winning_match[2]:
                winning_match = match         
        return winning_match
    else:
        return (None, None, None)
    
def match_coronae(bifucs, top_threshold, neighborhood_threshold=1):
    '''
    Takes the individual bifurcation lines and returns an ordered dict of (rank, x-axis point, corona -> bool): [(row,col)]
    '''

    matched_lines = set()
    coronae = OrderedDict()
    for i, ((n,v), pts) in enumerate(bifucs.iteritems()):

        # check that the current item hasn't already been matched against
        if n in matched_lines:
            continue
        else:
            # if the current number is not matched, then add it now
            matched_lines.add(n)

        # get the starting coordinates
        row_max, col_max = pts[0]
        row_min, col_min = pts[-1]

        # handle lines that traverse nearly the entire graph differently
        if row_max > top_threshold:
            coronae[(i, v, (pts[-1][1]))] = pts
            matched_lines.add(n)
            continue

        # This returns the first match against the current segment. The first should be the best.
        # pts_idx is the place to snip the tail off of the current line and glue the match and the current line together
        match_key, pts_idx, _ = get_best_match(matched_lines, pts, bifucs)

        if not match_key:
            coronae[(i, v, (pts[0][1]))] = pts
            continue
        match_pts = bifucs[match_key]

        # left and right here are just for clarity in logic, the parts may in fact be reversed in order (but it doesn't matter)
        left = pts[pts_idx:]
        right = match_pts

        # the 0-index has the max row for each part
        corona_max = left[0][0] if left[0][0] > right[0][0] else right[0][0]
        # now glue it all together
        corona = left[::-1] + right[:]

        coronae[(i, v, (corona[0][1],corona[-1][0]))] = corona

        # add the components to the matched set
        matched_lines.add(n)
        matched_lines.add(match_key[0])
        
    return coronae

In [7]:
# coronae = match_coronae(bifucations, 380, 1)
    
plt.figure(figsize=(14, 10))
# for n, (k, v) in enumerate(bifucations.items()):
#     rows, cols = zip(*v)
#     plt.plot(cols, rows, color='black')

for n, (k, v) in enumerate(bifucations.items()):
    rows, cols = zip(*v)
    plt.plot(cols, rows)
ax = plt.gca()
ax.set_ylim(ax.get_ylim()[::-1])
plt.setp(ax.get_xticklines(),visible=False)
plt.show()



In [5]:
import sys

In [6]:
sys.version_info


Out[6]:
sys.version_info(major=2, minor=7, micro=7, releaselevel='final', serial=0)

In [ ]: