def compute_mcm( g, estimate = None, pweight = 'weight' ):
    maxmean, arg_cycle = None, None
    for scc in nx.strongly_connected_component_subgraphs( g ):
        root = next( scc.nodes_iter() )
        scc_mcm, cycle = compute_mcm_component( scc, root, estimate, pweight )
        if scc_mcm is None:
            continue

        if maxmean is None or scc_mcm > maxmean:
            maxmean = scc_mcm
            arg_cycle = cycle

    forest = Forest()
    for scc in nx.strongly_connected_component_subgraphs( g, False ):
        if scc.number_of_edges() == 0:
            continue

        for ( v, w, scc_data ) in scc.edges_iter( data = True ):
            data = g.get_edge_data( v, w )

            # negate weight so that we can construct a longest paths tree for the current solution
            scc_data['w'] = maxmean - data.get( 'weight', 0 )

        root = w
        lpp_tree, _ = bfct.find_shortest_paths( scc, root, arg = 'w' )
        forest.add_forest( lpp_tree )

    return maxmean, arg_cycle, forest

def remove_all_cycles(g):
	n_scc = [x for x in nx.strongly_connected_component_subgraphs(g) if len(x)>1]
	if len(n_scc)>0:
		print "Found",len(n_scc),"cycles, largest one of size",len(n_scc[0])
	else:
		print "No cycles"
		return g
	# Color them 
	for ccc_idx  in range(len(n_scc)): 
		for node in n_scc[ccc_idx].nodes():
			g.node[node]['in_scc']=ccc_idx

	#sys.exit(0)
	g_p=g
	for i in range(1,1000): #Up to a thousand iteration

		g_n= remove_cycle(g_p)
		if g_n.number_of_edges() == g_p.number_of_edges():
			break
	#	assert("10_4" in g_n)
		g_p=g_n
		n_scc_i = [x for x in nx.strongly_connected_component_subgraphs(g_p) if len(x)>1]
		print "Still",len(n_scc_i),"cycles"

	g_trans=g_p
	# Assert graph is acyclic 
	assert(max([len(x) for x in nx.strongly_connected_components(g_trans)])==1)

	return g_trans

def cyclePlot(gexFile):
    DG = nx.DiGraph(nx.read_gexf(gexFile))
    
    #generate networkx friendly position format
    #dictionary keyed by node label with values being a float32 ndarray
    pos = dict()
    for i in range(1, len(DG.node)+1):
        xPos = DG.node[str(i)]['viz']['position']['x']
        yPos = DG.node[str(i)]['viz']['position']['y']
        pos[str(i)] = np.array([xPos,yPos])
    
    #nx.draw_networkx_edges(DG,pos,nodelist=DG.node.keys(),alpha=0.05,
    #                       arrows=True)
    nx.draw_networkx_nodes(DG,pos,nodelist=DG.node.keys(),
                       node_size=30,
                       node_color='grey',
                       alpha=0.4)
    nx.draw_networkx_edges(DG,pos,alpha=0.4,
                               arrows=True,
                               edge_color='k')
    plt.show()
    
    scc=nx.strongly_connected_component_subgraphs(DG)
    CG=scc[0];
    
    #show example
    nx.draw_networkx_nodes(CG,pos,nodelist=CG.node.keys(),
                       node_size=30,
                       node_color='c')
    nx.draw_networkx_edges(CG,pos,alpha=0.5,
                               arrows=True,
                               edge_color='r')

def main():

	
	parser = argparse.ArgumentParser()
	parser.add_argument('path', help = "target file or directory for summarization")
	parser.add_argument("--posseed", help="boosting seed for biased LexRank", default = None)
	parser.add_argument("--negseed", help="blocking seed for biased LexRank", default = None)
	parser.add_argument("--stopfile", help="file containing custom stopwords")
	parser.add_argument("-o", "--output", help = "output file name")
	parser.add_argument("-l", "--length", help = "summary length in lines", default = 10)
	parser.add_argument("--seed_posweight", help = "Weight for positive seed", default = 3)
	parser.add_argument("--seed_negweight", help = "Weight for negative seed", default = .0001)
	parser.add_argument("--ngrams", help = "N-gram number", default = 1)
	#normalization doesn't work due to being inherent within scoring method
	parser.add_argument("-n", "--is_norm", help = "Boolean flag for normalization", default = True)
	args = parser.parse_args()
	
	input_text = args.path
	pos_seed = args.posseed
	neg_seed = args.negseed
	stopfile = args.stopfile
	out_file = args.output
	sum_length = int(args.length)
	norm_flag = args.is_norm
	pos_weight = float(args.seed_posweight)
	neg_weight = float(args.seed_negweight)
	ngram = int(args.ngrams)
	corpus = Corpus(input_text).documents
	
	output_checker(out_file)

	if pos_seed == None and neg_seed == None:
		LR_method = 'unbiased'
		print LR_method
		[term_matrix, normalized] = TDM(corpus, pos_seed, neg_seed, stopfile, norm_flag, ngram).matrix
		pos_seed_vector = []
		neg_seed_vector = []
		
	else:
		LR_method = 'biased'
		if pos_seed == None:
			pos_seed = ''
		if neg_seed == None:
			neg_seed = ''
		
		[term_matrix, normalized, pos_seed_vector, neg_seed_vector] = TDM(corpus, pos_seed, neg_seed, stopfile, norm_flag, ngram).matrix
		corpus = corpus[2:]
		

	[scores,graph] = Graph(normalized, LR_method, pos_seed_vector, neg_seed_vector, pos_weight, neg_weight).sim_scores 
	#embed()
	largest = networkx.strongly_connected_component_subgraphs(graph)[0]
	A = networkx.to_agraph(largest)
	#A.node_attr.update(shape='point')
	A.node_attr.update(overlap='voronoi')
	A.layout(prog='sfdp')
	networkx.write_dot(largest, 'testgraph.gv')
	A.draw('butterflyv2.png')
	
	print_summary(corpus, scores, out_file, sum_length)

def strongly_connected_components():
    conn = sqlite3.connect("zhihu.db")     
    #following_data = pd.read_sql('select user_url, followee_url from Following where followee_url in (select user_url from User where agree_num > 50000) and user_url in (select user_url from User where agree_num > 50000)', conn)        
    following_data = pd.read_sql('select user_url, followee_url from Following where followee_url in (select user_url from User where agree_num > 10000) and user_url in (select user_url from User where agree_num > 10000)', conn)        
    conn.close()
    
    G = nx.DiGraph()
    cnt = 0
    for d in following_data.iterrows():
        G.add_edge(d[1][0],d[1][1])
        cnt += 1
    print 'links number:', cnt

    scompgraphs = nx.strongly_connected_component_subgraphs(G)
    scomponents = sorted(nx.strongly_connected_components(G), key=len, reverse=True)
    print 'components nodes distribution:', [len(c) for c in scomponents]
    
    #plot graph of component, calculate saverage_shortest_path_length of components who has over 1 nodes
    index = 0
    print 'average_shortest_path_length of components who has over 1 nodes:'
    for tempg in scompgraphs:
        index += 1
        if len(tempg.nodes()) != 1:
            print nx.average_shortest_path_length(tempg)
            print 'diameter', nx.diameter(tempg)
            print 'radius', nx.radius(tempg)
        pylab.figure(index)
        nx.draw_networkx(tempg)
        pylab.show()

    # Components-as-nodes Graph
    cG = nx.condensation(G)
    pylab.figure('Components-as-nodes Graph')
    nx.draw_networkx(cG)
    pylab.show()    

 def nontrivial_strongly_connected_components(self):
     ''' iterator over nontrivial strongly connected subgraphs\n
     nontrivial means it has more than one node or it is a node with a self loop '''
     scc = strongly_connected_component_subgraphs(self)
     for g in scc:
         if g.number_of_nodes()>1 or g.number_of_selfloops()>0:
             yield g

def test_dim_error():
    import sys
    authority_dict={}
    graph_file = '/home/michal/SALSA_files/tmp/real_run/middle_graph_authority'
    G_new = gm.read_graph_from_file(graph_file)
    isolates = nx.isolates(G_new)
    print 'num of isolates: '+str(len(isolates)); sys.stdout.flush()
    num_of_not_isolates = G_new.number_of_nodes() - len(isolates)
    authority_dict = {}
    classes = nx.strongly_connected_component_subgraphs(G_new)
    print 'num of classes including isolates: '+str(len(classes)); sys.stdout.flush()
    #remove classes of isolated nodes:   
    classes[:] = [ c for idx,c in enumerate(classes) if c.nodes()[0] not in isolates ]
    
    print 'num of classes NOT including isolates: '+str(len(classes)); sys.stdout.flush()
    for subG in classes:
        #print type(subG)
        out_file = ''.join(['/home/michal/SALSA_files/tmp/real_run/graph_',str(classes.index(subG))])
        gm.write_graph_to_file(subG, out_file)
        tmp_d = salsa.eig_calc(subG, normalize=num_of_not_isolates)    #power_iteration(subG)
    '''    
        for k,v in tmp_d.items():
            authority_dict[G.nodes()[k]] = v
        #print power_iteration(subG, tol=1.0e-10)
    for i in isolates:
        authority_dict[G.nodes()[i]] = 0
    #print authority_dict
    print '\n--- calc_salsa_per_class took: '+str(datetime.now()-startTime); sys.stdout.flush()'''
    return

def cycle_induced_subgraph(g):
    ''' Computes the subgraph that is maximally edge-induced by its cycles

    That is, every cycle that is in g is also in cycle_induced_subgraph(g).
    '''

    # gather all edges of SCCs
    return nx.DiGraph( ( edge for scc in nx.strongly_connected_component_subgraphs( g ) for edge in scc.edges_iter() ))

def remove_cycles(G):
    while not nx.is_directed_acyclic_graph(G):
        subgraphs = nx.strongly_connected_component_subgraphs(G)
        for subgraph in subgraphs:
            if subgraph.number_of_nodes() > 1:
                edge_index = random.randrange(subgraph.number_of_edges())
                edge = subgraph.edges()[edge_index]
                G.remove_edge(edge[0], edge[1])

def giant_component(net):
    if 'giant_component' in net.graph['already_computed'] :
        return net.graph['already_computed']['giant_component']
    else :
        if net.is_directed() :
            net.graph['already_computed']['giant_component'] = nx.strongly_connected_component_subgraphs(net)[0]
        else :
            net.graph['already_computed']['giant_component'] = nx.connected_component_subgraphs(net)[0]
        return net.graph['already_computed']['giant_component']

    def detect_deadlock(self):
        """
        Detects whether the system is in a deadlocked state, that is, is there a knot

        Note that this code is taken and adapted from the NetworkX Developer Zone Ticket #663 knot.py (09/06/2015)

            >>> from import_params import load_parameters
            >>> Q = Simulation(load_parameters('tests/datafortesting/logs_test_for_simulation/'))
            >>> nodes = ['A', 'B', 'C', 'D', 'E']
            >>> connections = [('A', 'D'), ('A', 'B'), ('B', 'E'), ('C', 'B'), ('E', 'C')]
            >>> for nd in nodes:
            ...     Q.digraph.add_node(nd)
            >>> for cnctn in connections:
            ...     Q.digraph.add_edge(cnctn[0], cnctn[1])
            >>> Q.detect_deadlock()
            True

            >>> Q = Simulation(load_parameters('tests/datafortesting/logs_test_for_simulation/'))
            >>> nodes = ['A', 'B', 'C', 'D']
            >>> connections = [('A', 'B'), ('A', 'C'), ('B', 'C'), ('B', 'D')]
            >>> for nd in nodes:
            ...     Q.digraph.add_node(nd)
            >>> for cnctn in connections:
            ...     Q.digraph.add_edge(cnctn[0], cnctn[1])
            >>> Q.detect_deadlock()
            False

            >>> Q = Simulation(load_parameters('tests/datafortesting/logs_test_for_simulation/'))
            >>> nodes = ['A', 'B']
            >>> for nd in nodes:
            ...     Q.digraph.add_node(nd)
            >>> Q.detect_deadlock()
            False
            >>> connections = [('A', 'A')]
            >>> for cnctn in connections:
            ...     Q.digraph.add_edge(cnctn[0], cnctn[1])
            >>> Q.detect_deadlock()
            True
        """
        knots = []
        for subgraph in nx.strongly_connected_component_subgraphs(self.digraph):
            nodes = set(subgraph.nodes())
            if len(nodes) == 1:
                n = nodes.pop()
                nodes.add(n)
                if set(self.digraph.successors(n)) == nodes:
                    knots.append(subgraph)
            else:
                for n in nodes:
                    successors = nx.descendants(self.digraph, n)
                    if not successors <= nodes:
                        break
                else:
                    knots.append(subgraph)
        if len(knots) > 0:
            return True
        return False

def basic_graph_process(_graph):
    # generate SCCs
    sccs = networkx.strongly_connected_component_subgraphs(_graph, copy=True)
    sccs = sorted(sccs, key=len, reverse=True)
    _graph.graph["sccs"] = sccs

    not_dag_edges = []

    for scc in sccs:
        not_dag_edges = not_dag_edges + scc.edges()

    dag_edges = [edge for edge in _graph.edges() if edge not in not_dag_edges]

    # give SCC index to each node, usefull ?????
    for i in range(0, len(sccs)):
        for node in sccs[i].nodes():
            _graph.node[node]["scc_index"] = i

    # find boundary nodes in each SCC
    for scc in sccs:
        scc.graph["inward_nodes"] = set()
        scc.graph["outward_nodes"] = set()

    # add inward and outward nodes
    for edge in dag_edges:
        scc_index = _graph.node[edge[0]]["scc_index"]
        scc = sccs[scc_index]
        scc.node[edge[0]].setdefault("outward_edges", [])
        scc.node[edge[0]]["outward_edges"].append(edge)
        scc.graph["outward_nodes"].add(edge[0])

        scc_index = _graph.node[edge[1]]["scc_index"]
        scc = sccs[scc_index]
        scc.node[edge[1]].setdefault("inward_edges", [])
        scc.node[edge[1]]["inward_edges"].append(edge)
        scc.graph["inward_nodes"].add(edge[1])

    # initial must be add to inward nodes, add final states as outward nodes
    initial_index = _graph.graph["initial"]
    sccs[_graph.node[initial_index]["scc_index"]].graph["inward_nodes"].add(initial_index)
    sccs[_graph.node[initial_index]["scc_index"]].node[initial_index].setdefault("inward_edges", [])

    final_indexes = _graph.graph["finals"]
    for final_index in final_indexes:
        scc_index = _graph.node[final_index]["scc_index"]
        scc = sccs[scc_index]
        scc.graph["outward_nodes"].add(final_index)
        scc.node[final_index].setdefault("outward_edges", [])


    # it's easier to operate on list ???
    for scc in sccs:
        scc.graph["inward_nodes"] = list(scc.graph["inward_nodes"])
        scc.graph["outward_nodes"] = list(scc.graph["outward_nodes"])

    return (sccs, dag_edges)

def get_major_strongly_component(graph):

    major_strongly_component_dict = {'amount_nodes': 0, 'subgraph':None}
    for component in nx.strongly_connected_component_subgraphs(graph):
        amount_nodes = len(component.nodes())
        if amount_nodes > major_strongly_component_dict['amount_nodes']:
            major_strongly_component_dict['amount_nodes'] = amount_nodes
            major_strongly_component_dict['subgraph'] = component

    return major_strongly_component_dict['subgraph']

    def __classifyEdges(self, tempG):
        """link classification
        :returns: @todo

        """
        Gd1, Gd2 = self.__turnMatchedGp2Gd(tempG)

        #get all the free nodes in Gp (self.BiG)
        freeNodes = filter(lambda res: res[1]['label']=='free', tempG.nodes(data=True))
        for node in freeNodes:
            nodeName = node[0]
            for Gd in [Gd1,Gd2]:
                flag, pathlist = self.__BFS_classifyEdges(Gd, nodeName)
                if flag:
                    for path in pathlist:
                        for i in range(len(path)-1):
                            Gd.edge[path[i]][path[i+1]]['mark'] = 'used'

        sccs = networkx.strongly_connected_component_subgraphs(Gd1)
        for subgraph in sccs:
            for link in subgraph.edges():
                Gd1.edge[link[0]][link[1]]['mark'] = 'used'
                Gd2.edge[link[1]][link[0]]['mark'] = 'used'

        for link in tempG.edges():
            flag0 = tempG.edge[link[0]][link[1]]['label'] == "matching"

            if link[0].endswith("\t+"):
                source = link[0]
                target = link[1]
            else:
                source = link[1]
                target = link[0]

            if flag0:
                flag1 = Gd1[source][target]['mark']=='unused'
                flag2 = Gd2[target][source]['mark']=='unused'
            else:
                flag1 = Gd1[target][source]['mark']=='unused'
                flag2 = Gd2[source][target]['mark']=='unused'

            # for random network
            source = source.strip("\t+")
            target = target.strip("\t-")

            if flag1 and flag2:
                if flag0:
                    self.DG.edge[source][target]['class'] = 'critical'
                else:
                    self.DG.edge[source][target]['class'] = 'redundant'
            else:
                self.DG.edge[source][target]['class'] = 'intermittent'

        pass

def myavgpathlength(G):
    try:
        apl =  nx.average_shortest_path_length(G)
        return [apl]
    except nx.NetworkXError as e:
        #this means graph is not connected
        if isinstance(G,nx.DiGraph):
		    return [nx.average_shortest_path_length(nx.strongly_connected_component_subgraphs(G)[0])]
        else:
            return [nx.average_shortest_path_length(nx.connected_component_subgraphs(G)[0])]
    except ZeroDivisionError as e:
        return [1]     

def linearize(filename):
    graph_name = filename.split('.')[0]+'.graphml'
    g = nx.read_graphml(graph_name)
    
    print nx.info(g)
    
    # get first strong connected component
    
    con = list(nx.strongly_connected_component_subgraphs(g))
    
    con.sort(key = lambda x:len(x), reverse = True)
    print [len(item) for item in con]

def diGraph_netX_stats(tmpDig):
	realSccs = 0
	scc = nx.strongly_connected_components(tmpDig)
	sccsg = nx.strongly_connected_component_subgraphs(tmpDig)
	actualScc = []
	for i in scc: 
		if len(i) > 1: actualScc.append(i)
	sccN= len(actualScc)
	selfLoops = tmpDig.number_of_selfloops()
	selfLoopsEgdes = tmpDig.selfloop_edges()
	realSccs = selfLoops + sccN 
	return actualScc, sccN, selfLoops, selfLoopsEgdes, realSccs, sccsg

	def filter_graph_for_strongly_connected_components(self, min_nodes=2): 
		"""
		Get strongly connected components in graph. 
		min_nodes : int 
			Return only connected components with a minimal number of min_nodes
		"""
		edges = [] 
		for g in nx.strongly_connected_component_subgraphs(self.graph): 
			if len(g.nodes(data=True)) >= min_nodes: 
				for e in g.edges(data=True): 
					edges.append(e)

		self.graph = nx.DiGraph(edges)

    def __init__(self, point_cloud, lines):
        # Simulates random order of lines
        random.shuffle(point_cloud)
        self.lines = lines
        self.point_cloud = point_cloud

        before = len(point_cloud)
        point_cloud = list(set(point_cloud))
        after = len(point_cloud)

        #assert before == after
        #assert len(point_cloud)/2 == len(lines)

        line_end_points = np.array(point_cloud)
        num_lines = len(line_end_points) / 2

        Y = cdist(np.array(line_end_points), np.array(line_end_points), 'euclidean')
        plt.figure(3)
        c = plt.imshow(Y, interpolation='none')
        plt.colorbar(c)
        plt.show()

        new = np.zeros(Y.shape)

        for i in range(len(Y)):
            row = Y[i]
            h = enumerate(row)
            h = sorted(h, key=lambda x: x[1])
            idx = [e[0] for e in h[:3]]
            for idx_i in idx:
                new[i][idx_i] = 1


        self.groups = []

        plt.figure(4)
        G = nx.from_numpy_matrix(new, create_using=nx.DiGraph())
        cols = itertools.cycle(["ro", "yo", "go", "mo", "bo", "co"])
        for subg in nx.strongly_connected_component_subgraphs(G):
            v = np.array([line_end_points[e] for e in subg.nodes()])
            print "Group", v.mean(), v.std()

            self.groups.append(v)


            c = cols.next()
            for element in subg.nodes():
                a, b = line_end_points[element]
                plt.plot(a, b, c)
        plt.show()

def calc_salsa_per_class_OLD(G, rank_type):
    '''
    # G is the original graph
    # rank_type = 'authority' or 'hub' 
    '''
    import scipy.sparse
    print '\t\t~~~~~~ calc_salsa_per_class ~~~~~~'; startTime = datetime.now(); sys.stdout.flush()

    # Authority:
    A=get_matrix(G,mat_type=rank_type,sparse=True)
    G_new = nx.DiGraph(A)
    '''DEBUG:
    out_file = ''.join(['/home/michal/SALSA_files/tmp/real_run/middle_graph_',str(rank_type)])
    gm.write_graph_to_file(G_new, out_file)'''
    #x=scipy.ones((n,1))/n  # initial guess
    isolates = nx.isolates(G_new)   # isolate node is a node with in_deg=out_deg=0
    print '--- calc_salsa_per_class: num of isolates- ',len(isolates),', out of- ',G_new.number_of_nodes(),' nodes (',float(len(isolates))/G_new.number_of_nodes(),'%)'; sys.stdout.flush()
    num_of_not_isolates = G_new.number_of_nodes() - len(isolates)
    scores_dict = {}
    tmpTime = datetime.now()
    classes = nx.strongly_connected_component_subgraphs(G_new)
    #print '--- calc_salsa_per_class: separate to classes took- '+str(datetime.now()-tmpTime); sys.stdout.flush(); tmpTime = datetime.now()
    #remove classes of isolated nodes:   
    classes[:] = [ c for idx,c in enumerate(classes) if c.nodes()[0] not in isolates ]
    #print '--- calc_salsa_per_class: clean classes from isolates took- ',datetime.now()-tmpTime; sys.stdout.flush(); 
    
    num_of_classes = 0
    domain_class_dict = {}
    for subG in classes:
        num_of_classes += 1
        '''DEBUG: 
        out_file = ''.join(['/home/michal/SALSA_files/tmp/real_run/graph_',str(classes.index(subG))])
        gm.write_graph_to_file(subG, out_file)'''
        tmp_d = eig_calc(subG, normalize=num_of_not_isolates)   
        #tmp_d = power_iteration(subG,max_iter=100,tol=1.0e-8,normalize=num_of_not_isolates,nstart=None)
        #tmp_d = power_iteration(subG, normalize=num_of_not_isolates, nstart=[v[G.n_attr.risk] for v in subG.nodes(data=True)])
        for k,v in tmp_d.items():
            d = G.nodes()[k]
            scores_dict[d] = v
            domain_class_dict[d] = num_of_classes
    print '--- calc_salsa_per_class: num of classes (NOT including isolates)- ',num_of_classes
    for i in isolates:
        d = G.nodes()[i]
        scores_dict[d] = 0
        domain_class_dict[d] = 0 # class zero represents the isolates
    #print authority_dict
    print '--- calc_salsa_per_class took: ',datetime.now()-startTime; sys.stdout.flush()

    return scores_dict, domain_class_dict