text:
import pandas as pd
import numpy as np
import seaborn as sns
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import copy
import numpy as np

def enriched_features(adata,f1,f2,fdr=None,x_or=5,x_p=20,**kwargs):

    pref='enrich_'+f1+'_vs_'+f2

    for item in ['.oddsratios','.p_adj.negLog10.signed']:
        if pref+item not in adata.uns:
            print("ERROR: result adata.uns['"+pref+item+"'] not found! \nPlease first run perturb.tl.enrich_features(adata,f1='"+f1+"',f2='"+f2+"')")
            return

    import seaborn as sns
    item='.oddsratios'
    g=sns.clustermap(adata.uns[pref+item],
                         vmin=-x_or,vmax=x_or,**kwargs)
    plt.title('log2 odds ratio\n('+f1+' vs '+f2+')',loc='center')
    ax = g.ax_heatmap
    ax.set_ylabel(f1)
    ax.set_xlabel(f2)

    item='.p_adj.negLog10.signed'
    g=sns.clustermap(adata.uns[pref+item],
                         vmin=-x_p,vmax=x_p,**kwargs)
    plt.title('signed -log10(p adj)\n('+f1+' vs '+f2+')',loc='center')
    ax = g.ax_heatmap
    ax.set_ylabel(f1)
    ax.set_xlabel(f2)

def perturbations_per_cell(adata_here,vmin=0,vmax=5,ax=None,pref='perturb',level='guide',copy=False):
    if ax==None:
        fig,plots=plt.subplots(1)
        ax=plots
    import copy
    vals=copy.deepcopy(adata_here.obs[pref+'.'+level+'.perturbations_per_cell'])
    vals[vals>vmax]=vmax
    vals[vals<vmin]=vmin

    ax.hist(vals,bins=range(int(vmin),int(vmax+1)),
                color='black')
    ax.set_xticks([vmin,vmax])
    ax.set_ylabel('Number of cells')
    ax.set_xlabel('Distinct '+level+'s\nper cell')
    ax.grid(False)
    plt.show()

def cells_per_perturbation(adata_here,ax=None,pref='perturb',level='guide',vmin=0,vmax=None):
    if ax==None:
        fig,plots=plt.subplots(1)
        ax=plots
    import copy
    for var in [pref+'.'+level+'.cells_per_perturbation',pref+'.'+level+'.cells_per_perturbation.singly_infected']:
        color='lightgray'
        label='all'
        if var==pref+'.'+level+'.cells_per_perturbation.singly_infected':
            color='black'
            label='single inf.'
        n, bins, patches=ax.hist(adata_here.uns[var]['Number of cells'],
                             color=color,
          cumulative=True,fill=False,histtype='step',label=label)
        patches[0].set_xy(patches[0].get_xy()[:-1])
    ax.grid(False)
    if vmin!=None and vmax!=None:
        ax.set_xlim(vmin,vmax)

    ax.set_xlabel('Cells per\n'+level)
    ax.set_ylabel('Number of\n'+level+'s')
    ax.legend(loc='lower right')

def ranked_cells_per_perturbation(adata_here,figwidth=2,figheight=10,pref='perturb',level='guide',vmin=0,vmax=None):
    if vmin!=None and vmax!=None:
        fig,plots=plt.subplots(1,2)
    else:
        fig,plots=plt.subplots(1,2,sharex=True)

    fig.set_size_inches(figwidth,figheight)

    import copy
    names=[pref+'.'+level+'.cells_per_perturbation.singly_infected',pref+'.'+level+'.cells_per_perturbation']
    for i in range(len(names)):
        df=adata_here.uns[names[i]]['Number of cells']

        color='lightgray'
        label='all'
        if names[i]==pref+'.'+level+'.cells_per_perturbation.singly_infected':
            color='black'
            label='single inf.'
            sorted_order=df.sort_values().index
        df=df.loc[sorted_order]

        plots[i].barh(df.index,df,color=color,label=label)
        plots[i].set_ylim(-0.5,df.shape[0]-0.5)
        plots[i].grid(False)
        if names[i]==pref+'.'+level+'.cells_per_perturbation.singly_infected':
            plots[i].set_yticklabels(sorted_order,fontsize=figheight/df.shape[0]*60)
        else:
            plots[i].set_yticklabels([])
        plots[i].set_title(label,fontsize=figwidth*5)
        plots[i].set_xlabel('Number\nof cells',fontsize=figwidth*5)
        if vmin!=None and vmax!=None:
            plots[i].set_xlim(vmin,vmax)

max_stars_repo_name: oursu/perturbseq

max_issues_repo_name: oursu/perturbseq

max_forks_repo_name: oursu/perturbseq


text:
"""
AUTHORS: Christopher J Burke
VERSION: 0.1
LICENSE:
MIT License

Copyright (c) 2021 Christopher Burke

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
"""

import gdb.types
from re import search as resrch
import numpy as np

# pyXt2Np is a subclass of gdb.Function class
class pyXt2Np (gdb.Function):
    """
    GDB convenience function that prints the dimensions and shape of an xtensor container using python. It also makes a numpy copy of the xtensor container as a global variable with the same variable name in the gdb interactive space, so that you can examine the data using python
    Usage: (gdb) p $xt2np("arr")
	where arr is an xt::xarray or xtensor or xtensor_fixed
	You must give the variable name as a string (e.g., surrounded in quotes)
    """
    def __init__ (self):
	# Call the __init__ of the parent class gdb.Function
    	# Which registers the convenience function name
    	super (pyXt2Np, self).__init__ ("xt2np")

    def invoke (self, var_name):
    	# Implementation of the xt2np convenience function within gdb
	# The first argument needs to be a string containing the variable
	#  name rather than the variable itself because
	#  we need the variable name later on to do method calls, save it as python global variable...
	# Name comes in with literal quotes so strip them away
	var_name_str = str(var_name).strip('"')
	# Get the actual gdb.Value for this variable
	var = gdb.parse_and_eval(var_name_str)

	# Determine if the variable passed is a valid xtensor type
	# Get the type of variable
	var_type = gdb.types.get_basic_type(var.type).name
	#print var_type
	# Look for xt::.*_container
	want_type = "xt::.*_container"
	if resrch(want_type, var_type):
	    # There was an xtensor container sent
	    # Confirm that the m_shape and m_storage fields are present
	    if gdb.types.has_field(var.type, 'm_shape') and gdb.types.has_field(var.type, 'm_storage'):
		# xtensor_fixed needs different treatment than xarray and xtensor to get shape and length
		want_type = "xt::xfixed_container"
		if not resrch(want_type, var_type):
		    # Handle cases of xarray and xtensor
                    # Get dimension from method call
		    dim = int(gdb.parse_and_eval('{}.dimension()'.format(var_name_str)))
		    shape_list = []
		    data_len = 1
		    for i in range(dim):
		    	curdim = int(gdb.parse_and_eval('{0}.shape()[{1:d}]'.format(var_name_str,i)).referenced_value())
		    	shape_list.append(curdim)
		    	data_len = data_len * curdim
		    # Get all the data from a gdb 'artifical array'
		    arr_data = gdb.parse_and_eval('*({0}.m_storage.p_begin)@{1:d}'.format(var_name_str, data_len))

		else:
		    # Handle case of xtensor_fixed
		    # The dimension shapes are constants and given in string representation of template
		    #  gdb.type.template_argument() can return those template arguments
                    shape_str = str(var.type.template_argument(1))
                    try:
		        found = resrch('<(.+?)>', shape_str).group(1)
		    except:
			return "Could not determine shape and size of a xtensor_fixed object"
		    shape_list = [int(x) for x in found.split(',')]
		    dim = len(shape_list)
		    data_len = 1
		    for i in range(dim):
			data_len = data_len * shape_list[i]
		    # Data is already in std::array
		    arr_data = gdb.parse_and_eval('{0}.m_storage._M_elems'.format(var_name_str))
		# python dictionary of accepted types go from C -> numpy dtypes
		nptypes = {'double':'d', 'float':'f', 'unsigned long':'L', 'unsigned int':'I', 'unsigned short':'H', 'long':'l','int':'i','short':'h'}
		gdb_type = arr_data[0].type.name
		if not gdb_type in nptypes:
			return "Could not find {} type equivalent to numpy".format(gdb_type)
		npUseType = nptypes[gdb_type]
		# make the numpy variable storage space
		np_arr_data = np.zeros((data_len,), dtype=npUseType)
		# iterate through the data array to load it into numpy array
		for i in range(data_len):
			np_arr_data[i] = arr_data[i]
		# Now reshape
		np_arr_fin = np_arr_data.reshape(shape_list)
		# Now make available global under same name as under gdb
		# According to google making global variable using exec and
	        #  making global variables in general is bad practice,
	        #  but here it is nonetheless
		com_str = "globals()['{}'] = np.copy(np_arr_fin)".format(var_name_str)
		exec(com_str)
		# Use numpy's builin pretty printing to show the contents of xtensor
		exec("print({})".format(var_name_str))
		# Return the make stats about xtensor
	        return "xtensor {0} Dimen: {1:d} Shape: {2} Len: {3:d} GDBType: {4} NPType: {5}-{6}".format(var_name_str,dim, str(shape_list), data_len, gdb_type, npUseType, np_arr_fin.dtype)
	else:
	    # There was not an xtensor container as argument return
	    return "The argument to $xt2np is not an xtensor container"

# Need to instantiate the class
pyXt2Np()

ext: py

reason: import numpy


text:
import sys
import nltk
nltk.download(['punkt', 'wordnet', 'stopwords'])
from nltk.tokenize import word_tokenize
from nltk.corpus import stopwords
from nltk.stem import WordNetLemmatizer

import re
import numpy as np
import pandas as pd
import pickle
import sklearn
from sqlalchemy import create_engine

from sklearn.metrics import classification_report, make_scorer
from sklearn.model_selection import GridSearchCV
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import train_test_split
from sklearn.pipeline import Pipeline
from sklearn.feature_extraction.text import CountVectorizer, TfidfTransformer
from sklearn.multioutput import MultiOutputClassifier
from sklearn.multiclass import OneVsRestClassifier
from sklearn.svm import LinearSVC

def load_data(database_filepath):

    '''
    FUNCTION:
        Load database file containing the cleaned data

    INPUTS:
        database_filepath - database file containing the cleaned data

    OUTPUTS
        X - The messages from the database
        y - the categories the messages fall in marked 0 or 1
        y.columns - the names of the categories

    '''

    engine = create_engine('sqlite:///{}'.format(database_filepath))
    df = pd.read_sql_table(database_filepath, engine)
    X = df['message']
    y = df.iloc[:,4:]
    return X, y, y.columns

def tokenize(text):

    '''
    FUNCTION:
        clean the input text, simplify it, and tokenize it

    INPUTS:
        text - a string to clean and tokenize

    OUTPUTS
        clean_tokens - a list of tokenized words from text
    '''

    #remove punctuation and convert to lowercase
    message = re.sub(r'[^a-zA-Z0-9]', ' ', text.lower())

    #tokenize
    tokens = word_tokenize(message)

    #lemmatization
    lemmatizer = WordNetLemmatizer()

    clean_tokens = []
    for tok in tokens:
        if tok not in stopwords.words('english'):
            clean_tok = lemmatizer.lemmatize(tok).strip()
            clean_tokens.append(clean_tok)

    return clean_tokens

def gscv_scores(y_test, y_pred):

    '''
    Function: calculate the f-1 score average of all classes to use in
              GridSearch evaluation
    Inputs:
        -y_test: outputs of test data set from train_test_split()
        -y_pred: pipeline predicted output based on predicted x from
                 train_test_split()
    Outputs:
        -average f-1 score of all classes to determine if a parameter improved
         predictibility
    '''
    y_pred = pd.DataFrame(y_pred, columns = y_test.columns)
    report = pd.DataFrame()

    for col in y_test.columns:
        class_rep = classification_report(y_test.loc[:,col], y_pred.loc[:,col],
                                          output_dict=True, zero_division=0)
        scoring_df = pd.DataFrame.from_dict(class_rep)

        #print(scoring_df)

        #drop un-needed columns and rows
        scoring_df = scoring_df[['macro avg']]
        scoring_df.drop(index='support', inplace=True)
        scoring_df = scoring_df.transpose()
        report = report.append(scoring_df, ignore_index=True)

    report.index = y_test.columns
    return report['f1-score'].mean()

def build_model():

    '''
    FUNCTION: create and return a pipeline
    '''
    pipeline = Pipeline([('vect', CountVectorizer(tokenizer=tokenize)),
                     ('tfidf', TfidfTransformer()),
                     ('clf', MultiOutputClassifier(RandomForestClassifier()))])

    #comment different parameters to decrease run time
    parameters = {
        #'vect__ngram_range': ((1,1), (1,2)),
        #'vect__max_df': (0.5, 0.75, 1.0),
        'vect__max_features': (None, 5000, 10000),
        #'tfidf__use_idf': (True, False),
        #'tfidf__norm': ('l1', 'l2'),
        'clf__estimator__n_estimators': [10, 25, 50],
        #'clf__estimator__min_samples_split': [2, 4, 6],
        'clf__estimator__bootstrap': (True, False)
    }

    scorer = make_scorer(gscv_scores)

    cv = GridSearchCV(pipeline, param_grid=parameters, scoring=scorer, cv=3, verbose=4)
    return cv

def evaluate_model(model, X_test, Y_test, category_names):

    '''
    FUNCTION:
        Evaluate the effectiveness of the pipeline in modeling the data

    INPUTS:
        model - pipeline fit to the data
        X_test - test set of messages from the train_test_split
        Y_test - test set of categories from the train_test_split
        category_names - names of the categories the messages can fit in

    OUTPUTS
        report - the precision, recall, and f-1 scores for each category
        scores - the average precision, recal, and f-1 scores for the data set
    '''

    y_pred = model.predict(X_test)
    y_pred = pd.DataFrame(y_pred, columns=category_names)

    report = pd.DataFrame()

    for col in Y_test.columns:

        class_rep = classification_report(Y_test.loc[:,col],
                                          y_pred.loc[:,col],
                                          output_dict=True,
                                          zero_division=0)
        scoring_df = pd.DataFrame.from_dict(class_rep)

        scoring_df = scoring_df[['macro avg']]
        scoring_df.drop(index='support', inplace=True)
        scoring_df = scoring_df.transpose()
        report = report.append(scoring_df, ignore_index=True)

    report.index = category_names
    scores = {'mean precision score':report['precision'].mean(),
              'mean recall score':report['recall'].mean(),
              'mean f-1 score':report['f1-score'].mean()}

    print(report)
    print(scores)

def save_model(model, model_filepath):

    '''
    FUNCTION:
        save the model as a pickle file

    INPUTS:
        model - the pipeline trained to the data
        model_filepath - where to save the file

    OUTPUTS
        None
    '''

    pickle.dump(model, open(model_filepath, 'wb'))

def main():
    if len(sys.argv) == 3:
        database_filepath, model_filepath = sys.argv[1:]
        print('Loading data...\n    DATABASE: {}'.format(database_filepath))
        X, Y, category_names = load_data(database_filepath)
        X_train, X_test, Y_train, Y_test = train_test_split(X, Y, test_size=0.2)

        print('Building model...')
        model = build_model()

        print('Training model...')
        model.fit(X_train, Y_train)

        print('Evaluating model...')
        evaluate_model(model, X_test, Y_test, category_names)

        print('Saving model...\n    MODEL: {}'.format(model_filepath))
        save_model(model, model_filepath)

        print('Trained model saved!')

    else:
        print('Please provide the filepath of the disaster messages database '\
              'as the first argument and the filepath of the pickle file to '\
              'save the model to as the second argument. \n\nExample: python '\
              'train_classifier.py ../data/DisasterResponse.db classifier.pkl')

if __name__ == '__main__':
    main()

max_stars_repo_name: jadefreese/DisasterResponse

max_stars_repo_licenses: FTL CNRI-Python blessing

max_issues_repo_name: jadefreese/DisasterResponse

max_issues_repo_licenses: FTL CNRI-Python blessing

max_forks_repo_name: jadefreese/DisasterResponse

max_forks_repo_licenses: FTL CNRI-Python blessing


text:
from torch.autograd import Variable
from tqdm import tqdm
from torchvision import datasets, models, transforms
import torch
import torch.nn as nn
import torch.optim as optim
from torch.optim import lr_scheduler
import numpy as np
import time
import os
import argparse
import torchvision
import copy
import sys
sys.path.append("../..")
from utils.data_processing import get_data_transforms, load_data, feature_eval_prep
from modules.learning import train_model
from modules.model_config import model_take_lower_layers
from modules.models.Inception_V3_trimmed import inception_v3_trimmed
from torchvision.models import inception_v3

def node_training(data_dir='./data/', restore=1, model_path='./checkpoint_i/base_model', freeze_layers=1, n_class=2,
                  batch_size=12, epochs = 25, device_ids=[0,1], output_dir = './checkpoint_o/', new_model_name='new_model.t7' ):
	### This script is to train model on the raw image data from clients. The dataset
	### should be arranged in train and val folders before running this code.

	## Loading the dataloaders -- Make sure that the data is saved in following way
	"""
	data/
	  - train/
	      - class_1 folder/
		  - img1.png
		  - img2.png
	      - class_2 folder/
	      .....
	      - class_n folder/
	  - val/
	      - class_1 folder/
	      - class_2 folder/
	      ......
	      - class_n folder/
	"""
	#parser = argparse.ArgumentParser(description='PyTorch inception Training')
	#parser.add_argument('--img_dir', default='./refold1/', type=str, help='')
	#parser.add_argument('--output_dir', default='./checkpoint/', type=str, help='')
	#parser.add_argument('--restore', default=1, type=int, help='whether restore the stored model or not')
	#parser.add_argument('--final_output', default=2, type=int, help='the number of the output of the final layer')
	#parser.add_argument('--batch_size', default=12, type=int)
	#parser.add_argument('--epochs', default=30, type=int)

	#args = parser.parse_args()
	#data_dir = args.img_dir
	#restore=args.restore
	#n_class=args.final_output
	#batch_size=args.batch_size
	#epochs=args.epochs
	#output_dir=args.output_dir

	data_dir = data_dir
	restore = restore
	model_path = model_path
	freeze_layers=freeze_layers
	n_class = n_class
	batch_size = batch_size
	epochs = epochs
	output_dir=output_dir
	new_model_name = new_model_name

	input_shape = 299
	mean = [0.485, 0.456, 0.406]
	std = [0.229, 0.224, 0.225]
	scale = 299
	input_shape = 299
	use_parallel = True
	use_gpu = True

	data_transforms = get_data_transforms(mode="train", model_origin="Inception V3",
		                  mean=None, std=None, scale=None, input_shape=input_shape)

	dataloaders, dataset_sizes, class_names, class_2_idx = load_data(data_dir,
		                                                         data_transforms,
		                                                         mode = 'train',
		                                                         batch_size=batch_size,
		                                                         num_workers=None)

	if(restore):
	    print('==> Resuming from checkpoint..')
	    assert os.path.isdir('checkpoint_i'), 'Error: no checkpoint directory found!'
	    model_conv=torch.load(model_path)
	    model_conv=model_conv.module

	else:
	    model_conv = inception_v3(pretrained=True)

	"""
	freeze_layers=1 is a restricted form of transfer learning with only the fully connected layers open for update. If this 	variable is 0, then all parameters of the model is updated in the training.
	"""
	if freeze_layers:
	    for i, param in model_conv.named_parameters():
	        param.requires_grad = False

	# Since imagenet as 1000 classes , We need to change our last layer according to the number of classes we have,
	num_ftrs = model_conv.fc.in_features
	model_conv.fc = nn.Linear(num_ftrs, n_class)

	# Stage-2 , Freeze all the layers till "Conv2d_4a_3*3"
	ct = []
	for name, child in model_conv.named_children():
	    # Mixed_7c and higher are the layers to be retrained
	    if "Mixed_7c" in ct:
	    	for params in child.parameters():
	    		params.requires_grad = True
	    ct.append(name)

	#if use_gpu:
	#    model_conv.cuda()

	#if use_parallel:
	#    print("[Using all the available cores]")
	#    model_conv = nn.DataParallel(model_conv, device_ids=[0, 1])

	if use_gpu:
	    model_conv.cuda()
	    print("[Using all the available GPUs]")
	    if use_parallel:
	    	model_conv = nn.DataParallel(model_conv, device_ids=device_ids)
	    	model_conv.cuda()
	    print("[Using CrossEntropyLoss...]")
	    criterion = nn.CrossEntropyLoss().cuda()
	else:
	    if not use_parallel:
	    	model_conv = nn.DataParallel(model_conv, device_ids=device_ids)
	    print("[Using CrossEntropyLoss...]")
	    criterion = nn.CrossEntropyLoss()

	print("[Using small learning rate with momentum...]")
	optimizer = optim.SGD(list(filter(lambda p: p.requires_grad, model_conv.parameters())), lr=0.005, momentum=0.9)

	print("[Creating Learning rate scheduler...]")
	scheduler = lr_scheduler.StepLR(optimizer, step_size=7, gamma=0.1)

	print("[Training the model begins ....]")
	train_model(model_conv, dataloaders, dataset_sizes, criterion, optimizer,
		    scheduler, use_gpu, epochs, output_dir, output_name=new_model_name, mixup = False, alpha = 0.1)

	# Model is saved in output_dir

	# Then send new model to the clients if necessary.

#if __name__=__main__:
#	node_training()


text:
# -*- coding: utf-8 -*-
"""
Created on Thu Jun 8

@author: rohankoodli
"""

import numpy as np
import os
from readData import read_movesets_pid, read_structure
from encodeRNA import encode_movesets, encode_structure
import tensorflow as tf
from tensorflow.contrib import rnn
import pickle
# enc0 = np.array([[[[1,2,3,4],[0,1,0,1],[-33,0,0,0]],[[1,2,3,4],[0,1,1,0],[-23,0,0,0]]],[[[3,3,3,3],[0,0,0,0],[2,0,0,0]],[[1,1,1,0],[1,0,1,0],[-23,0,0,0]]]])
# ms0 = np.array([[[2,1],[4,3]],[[1,6],[2,9]]])
# enc = np.array([[[1,2,3,4],[0,1,0,1],[1,1,1,1],[-3,0,0,0]],[[4,3,2,1],[1,0,1,0],[0,0,0,0],[9,0,0,0]]])
# out = np.array([[4,2],[3,3]])

real_X = pickle.load(open(os.getcwd()+'/pickles/X-6892348','rb'))
real_y = pickle.load(open(os.getcwd()+'/pickles/y-6892348','rb'))

real_X_9 = np.array(real_X[0:1500]).reshape([-1,340])
real_y_9 = np.array(real_y[0:1500])
test_real_X = np.array(real_X[1500:1600]).reshape([-1,340])
test_real_y = np.array(real_y[1500:1600])

#real_X_9, test_real_X, real_y_9, test_real_y = np.array(train_test_split(real_X[0:500],real_y[0:500],test_size=0.2))
#real_X_9, test_real_X, real_y_9, test_real_y = np.array(real_X_9).reshape([-1,340]), np.array(test_real_X).reshape([-1,340]), np.array(real_y_9), np.array(test_real_y)

# enc0 = np.array([[[1,2,3,4],[0,1,0,1],[-33,0,0,0],[1,1,1,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]],[[2,3,3,2],[0,0,0,0],[9,0,0,0],[0,0,0,1]]])
# ms0 = np.array([[1,6],[2,7],[2,7],[2,7],[2,7],[2,7],[2,7],[2,7],[2,7]])
# ms0 = np.array([[1,0,0,0],[0,1,0,0],[0,1,0,0],[0,1,0,0],[0,1,0,0],[0,1,0,0],[0,1,0,0],[0,1,0,0],[0,1,0,0]]) # just base
#
# test_enc0 = np.array([[[2,3,3,2],[0,0,0,0],[6,0,0,0],[0,0,1,1]],[[1,2,3,4],[0,1,0,1],[-33,0,0,0],[1,1,1,1]]])
# test_ms0 = np.array([[4,20],[3,15]])
# test_ms0 = np.array([[0,0,0,1],[1,0,0,0]]) # just base

num_epochs = 10
n_classes = 4
batch_size = 1500 # same as training size
chunk_size = 4
n_chunks = 85
rnn_size = 128
TRAIN_KEEP_PROB = 0.5

x = tf.placeholder('float',[None,n_chunks,chunk_size]) # 16 with enc0
y = tf.placeholder('float')
keep_prob = tf.placeholder('float')

# enc = enc0.reshape([-1,16])
# ms = ms0#.reshape([-1,4])
#
# test_enc = test_enc0.reshape([-1,16])
# test_ms = test_ms0
#e1 = tf.reshape(enc0,[])

def recurrentNeuralNet(x):
    layer = {'weights':tf.Variable(tf.random_normal([rnn_size,n_classes])),
            'biases':tf.Variable(tf.random_normal([n_classes]))}

    x = tf.transpose(x, [1,0,2])
    x = tf.reshape(x, [-1,chunk_size])
    x = tf.split(x, n_chunks, 0)

    lstm_cell = rnn.BasicLSTMCell(rnn_size)

    lstm_cell = rnn.DropoutWrapper(lstm_cell, output_keep_prob=keep_prob)

    outputs, states = rnn.static_rnn(lstm_cell, x, dtype=tf.float32)

    ol = tf.matmul(outputs[-1], layer['weights']) + layer['biases']

    return ol

def train(x):
    prediction = recurrentNeuralNet(x)
    #print prediction
    cost = tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=prediction,labels=y))
    optimizer = tf.train.AdamOptimizer().minimize(cost) # learning rate = 0.001

    # cycles of feed forward and backprop

    with tf.Session() as sess:
        sess.run(tf.global_variables_initializer())

        for epoch in range(num_epochs):
            epoch_loss = 0
            for _ in range(int(real_X_9.shape[0])/batch_size):#mnist.train.num_examples/batch_size)): # X.shape[0]
                randidx = np.random.choice(real_X_9.shape[0], batch_size, replace=False)
                epoch_x,epoch_y = real_X_9[randidx,:],real_y_9[randidx,:] #mnist.train.next_batch(batch_size) # X,y
                epoch_x = epoch_x.reshape((batch_size, n_chunks, chunk_size))
                _,c = sess.run([optimizer,cost],feed_dict={x:epoch_x,y:epoch_y, keep_prob: TRAIN_KEEP_PROB})
                epoch_loss += c
            print 'Epoch', epoch + 1, 'completed out of', num_epochs, '\nLoss:',epoch_loss,'\n'

        correct = tf.equal(tf.argmax(prediction,1),tf.argmax(y,1))
        accuracy = tf.reduce_mean(tf.cast(correct,'float'))

        print 'Accuracy', accuracy.eval(feed_dict={x:test_real_X.reshape((-1, n_chunks, chunk_size)), y:test_real_y, keep_prob: 1.0}) #X, y #mnist.test.images, mnist.test.labels

        # data for tensorboard
        writer = tf.summary.FileWriter(os.getcwd()+'/tensorboard/baseRNN-1500-100-128-10')
        writer.add_graph(sess.graph)
        '''
        Run this:
        tensorboard --logdir=tensorboard/baseRNN-SPECIFICATIONS --debug
        '''

train(x)

max_stars_repo_name: eternagame/EternaBrain

max_issues_repo_name: Eternagame/EternaBrain

max_forks_repo_name: Eternagame/EternaBrain


text:
# MIT License
# code by Soohwan Kim @sooftware

import wandb
import os
import torch
import logging
import numpy as np
import pandas as pd
from tqdm import tqdm

from utils import code_to_result

logger = logging.getLogger(__name__)

class Trainer:
    def __init__(
            self,
            model,
            train_loader,
            valid_loader,
            optimizer,
            scheduler,
            num_epochs,
            device,
            tokenizer,
            gradient_clip_val: float,
            accumulate_grad_batches: int,
            sep_token: str,
            log_every: int = 20,
            save_dir: str = 'ckpt',
    ) -> None:
        super(Trainer, self).__init__()
        self.model = model
        self.train_loader = train_loader
        self.valid_loader = valid_loader
        self.optimizer = optimizer
        self.scheduler = scheduler
        self.log_every = log_every
        self.num_epochs = num_epochs
        self.tokenizer = tokenizer
        self.eos_token = tokenizer.eos_token
        self.sep_token = sep_token
        self.device = device
        self.gradient_clip_val = gradient_clip_val
        self.accumulate_grad_batches = accumulate_grad_batches
        self.save_dir = save_dir

    def _get_lr(self):
        for param_group in self.optimizer.param_groups:
            return param_group['lr']

    def _save_result(self, inputs, codes, predicts, answers, epoch):
        pd.DataFrame({
            "inputs": inputs,
            "generate_codes": codes,
            "predicts": predicts,
            "answers": answers,
        }).to_csv(f'{epoch}_result.csv', index=False)

    def _train_epoch(self):
        total_loss = 0.0

        self.model.train()

        for step, (inputs, attention_masks, labels, answers) in enumerate(tqdm(self.train_loader)):
            inputs = inputs.to(self.device)
            attention_masks = attention_masks.to(self.device)
            labels = labels.to(self.device)

            batch_size = inputs.size(0)

            outputs = self.model(input_ids=inputs, attention_mask=attention_masks, labels=labels)
            loss = outputs.loss

            loss = loss / self.accumulate_grad_batches

            self.optimizer.zero_grad()
            loss.backward()
            total_loss += loss.item()

            if step % self.accumulate_grad_batches == 0:
                torch.nn.utils.clip_grad_norm_(self.model.parameters(), self.gradient_clip_val)
                self.optimizer.step()
                self.scheduler.step()

            if step % self.log_every == 0:
                train_mean_loss = total_loss / (step + 1)
                wandb.log({"lr": self._get_lr(),
                           "train_loss": train_mean_loss,
                           "batch_size": batch_size})

    def _get_accuracy(self, predicts, answers):
        accuracy = 0.0
        num = len(predicts)

        for pred, answer in zip(predicts, answers):
            if pred == answer:
                accuracy += 1

        accuracy /= num

        return accuracy

    def _validate(self, epoch):
        all_inputs = list()
        all_codes = list()
        all_predicts = list()
        all_answers = list()

        self.model.eval()

        with torch.no_grad():
            for step, (input_ids, attention_masks, labels, answers) in enumerate(tqdm(self.valid_loader)):
                input_ids = input_ids.to(self.device)
                attention_masks = attention_masks.to(self.device)

                generated = self.model.generate(
                    input_ids=input_ids,
                    attention_mask=attention_masks,
                    max_length=300,
                )
                generated = self.tokenizer.batch_decode(generated, skip_special_tokens=False)

                for g, a in zip(generated, answers):
                    splited = g.split(self.sep_token)
                    input_var, code = splited[0], splited[1].split(self.eos_token)[0]
                    code = code.replace('<pad>', '').replace('<s>', '').replace('<unk>', '').strip()
                    try:
                        predict = code_to_result(code)
                    except:
                        predict = None

                    all_inputs.append(input_var)
                    all_codes.append(code)
                    if predict is not None:
                        all_predicts.append(predict)
                    else:
                        all_predicts.append(np.nan)
                    all_answers.append(a)

        accuracy = self._get_accuracy(all_predicts, all_answers)
        wandb.log({"valid_accuracy": accuracy})
        self._save_result(all_inputs, all_codes, all_predicts, all_answers, epoch)

        if not os.path.exists(self.save_dir):
            os.mkdir(self.save_dir)
        self.model.save_pretrained(os.path.join(self.save_dir, f"epoch_{epoch}"))

    def fit(self):
        for epoch in range(self.num_epochs):
            logger.info(f"Epoch {epoch} start..")
            self._train_epoch()

            logger.info(f"Validate..")
            self._validate(epoch)

        self.model.save_pretrained(os.path.join(self.save_dir, "LAST_CHECKPOINT"))

max_stars_repo_name: sooftware/nlp-tasks

max_issues_repo_name: sooftware/nlp-tasks

max_forks_repo_name: sooftware/nlp-tasks


text:
# Copyright 2019 The Cirq Developers
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     https://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
"""Tests for measures."""
import numpy as np
import pytest

import cirq

N = 15
VEC1 = cirq.testing.random_superposition(N)
VEC2 = cirq.testing.random_superposition(N)
MAT1 = cirq.testing.random_density_matrix(N)
MAT2 = cirq.testing.random_density_matrix(N)
U = cirq.testing.random_unitary(N)

def test_fidelity_symmetric():
    np.testing.assert_allclose(cirq.fidelity(VEC1, VEC2), cirq.fidelity(VEC2, VEC1))
    np.testing.assert_allclose(cirq.fidelity(VEC1, MAT1), cirq.fidelity(MAT1, VEC1))
    np.testing.assert_allclose(
        cirq.fidelity(cirq.density_matrix(MAT1), MAT2),
        cirq.fidelity(cirq.density_matrix(MAT2), MAT1),
    )

def test_fidelity_between_zero_and_one():
    assert 0 <= cirq.fidelity(VEC1, VEC2) <= 1
    assert 0 <= cirq.fidelity(VEC1, MAT1) <= 1
    assert 0 <= cirq.fidelity(cirq.density_matrix(MAT1), MAT2) <= 1

def test_fidelity_invariant_under_unitary_transformation():
    np.testing.assert_allclose(
        cirq.fidelity(cirq.density_matrix(MAT1), MAT2),
        cirq.fidelity(cirq.density_matrix(U @ MAT1 @ U.T.conj()), U @ MAT2 @ U.T.conj()),
    )

def test_fidelity_commuting_matrices():
    d1 = np.random.uniform(size=N)
    d1 /= np.sum(d1)
    d2 = np.random.uniform(size=N)
    d2 /= np.sum(d2)
    mat1 = cirq.density_matrix(U @ np.diag(d1) @ U.T.conj())
    mat2 = U @ np.diag(d2) @ U.T.conj()

    np.testing.assert_allclose(
        cirq.fidelity(mat1, mat2, qid_shape=(15,)), np.sum(np.sqrt(d1 * d2)) ** 2
    )

def test_fidelity_known_values():
    vec1 = np.array([1, 1j, -1, -1j]) * 0.5
    vec2 = np.array([1, -1, 1, -1], dtype=np.complex128) * 0.5
    vec3 = np.array([1, 0, 0, 0], dtype=np.complex128)
    tensor1 = np.reshape(vec1, (2, 2))
    mat1 = cirq.density_matrix(np.outer(vec1, vec1.conj()))
    mat2 = cirq.density_matrix(np.outer(vec2, vec2.conj()))
    mat3 = 0.3 * mat1.density_matrix() + 0.7 * mat2.density_matrix()

    np.testing.assert_allclose(cirq.fidelity(vec1, vec1), 1)
    np.testing.assert_allclose(cirq.fidelity(vec2, vec2), 1)
    np.testing.assert_allclose(cirq.fidelity(vec1, vec3), 0.25)
    np.testing.assert_allclose(cirq.fidelity(vec1, tensor1), 1)
    np.testing.assert_allclose(cirq.fidelity(tensor1, vec1), 1)
    np.testing.assert_allclose(cirq.fidelity(mat1, mat1), 1)
    np.testing.assert_allclose(cirq.fidelity(mat2, mat2), 1)
    np.testing.assert_allclose(cirq.fidelity(vec1, mat1), 1)
    np.testing.assert_allclose(cirq.fidelity(mat2, vec2), 1)
    np.testing.assert_allclose(cirq.fidelity(vec1, vec2), 0)
    np.testing.assert_allclose(cirq.fidelity(vec1, mat2), 0)
    np.testing.assert_allclose(cirq.fidelity(mat1, vec2), 0)
    np.testing.assert_allclose(cirq.fidelity(vec1, mat3), 0.3)
    np.testing.assert_allclose(cirq.fidelity(tensor1, mat3), 0.3)
    np.testing.assert_allclose(cirq.fidelity(mat3, tensor1), 0.3)
    np.testing.assert_allclose(cirq.fidelity(mat3, vec2), 0.7)

def test_fidelity_numpy_arrays():
    vec1 = np.array([1, 0, 0, 0, 0, 0, 0, 0], dtype=np.complex64)
    vec2 = np.array([1, 0, 0, 0], dtype=np.complex64)
    tensor1 = np.reshape(vec1, (2, 2, 2))
    tensor2 = np.reshape(vec2, (2, 2))
    mat1 = np.outer(vec1, vec1.conj())

    np.testing.assert_allclose(cirq.fidelity(vec1, vec1), 1)
    np.testing.assert_allclose(cirq.fidelity(vec1, tensor1), 1)
    np.testing.assert_allclose(cirq.fidelity(tensor1, vec1), 1)
    np.testing.assert_allclose(cirq.fidelity(vec1, mat1), 1)
    np.testing.assert_allclose(cirq.fidelity(tensor1, mat1), 1)
    np.testing.assert_allclose(cirq.fidelity(tensor2, tensor2, qid_shape=(2, 2)), 1)
    np.testing.assert_allclose(cirq.fidelity(mat1, mat1, qid_shape=(8,)), 1)

    with pytest.raises(ValueError, match='dimension'):
        _ = cirq.fidelity(vec1, vec1, qid_shape=(4,))

    with pytest.raises(ValueError, match='Mismatched'):
        _ = cirq.fidelity(vec1, vec2)

    with pytest.raises(ValueError, match='ambiguous'):
        _ = cirq.fidelity(tensor2, tensor2)

    with pytest.raises(ValueError, match='ambiguous'):
        _ = cirq.fidelity(mat1, mat1)

def test_fidelity_ints():
    assert cirq.fidelity(3, 4) == 0.0
    assert cirq.fidelity(4, 4) == 1.0

    with pytest.raises(ValueError, match='non-negative'):
        _ = cirq.fidelity(-1, 2)
    with pytest.raises(ValueError, match='maximum'):
        _ = cirq.fidelity(4, 1, qid_shape=(2,))
    with pytest.raises(ValueError, match='maximum'):
        _ = cirq.fidelity(1, 4, qid_shape=(2,))

def test_fidelity_product_states():
    a, b = cirq.LineQubit.range(2)

    np.testing.assert_allclose(
        cirq.fidelity(cirq.KET_ZERO(a) * cirq.KET_ZERO(b), cirq.KET_ZERO(a) * cirq.KET_ZERO(b)), 1.0
    )
    np.testing.assert_allclose(
        cirq.fidelity(cirq.KET_ZERO(a) * cirq.KET_ZERO(b), cirq.KET_ZERO(a) * cirq.KET_ONE(b)),
        0.0,
        atol=1e-7,
    )
    np.testing.assert_allclose(
        cirq.fidelity(cirq.KET_ZERO(a) * cirq.KET_ZERO(b), cirq.KET_ZERO(a) * cirq.KET_PLUS(b)), 0.5
    )
    np.testing.assert_allclose(
        cirq.fidelity(cirq.KET_ONE(a) * cirq.KET_ONE(b), cirq.KET_MINUS(a) * cirq.KET_PLUS(b)), 0.25
    )
    np.testing.assert_allclose(
        cirq.fidelity(cirq.KET_MINUS(a) * cirq.KET_PLUS(b), cirq.KET_MINUS(a) * cirq.KET_PLUS(b)),
        1.0,
    )
    np.testing.assert_allclose(
        cirq.fidelity(cirq.KET_MINUS(a) * cirq.KET_PLUS(b), cirq.KET_PLUS(a) * cirq.KET_MINUS(b)),
        0.0,
        atol=1e-7,
    )

    with pytest.raises(ValueError, match='Mismatched'):
        _ = cirq.fidelity(cirq.KET_MINUS(a), cirq.KET_PLUS(a) * cirq.KET_MINUS(b))
    with pytest.raises(ValueError, match='qid shape'):
        _ = cirq.fidelity(
            cirq.KET_MINUS(a) * cirq.KET_PLUS(b),
            cirq.KET_PLUS(a) * cirq.KET_MINUS(b),
            qid_shape=(4,),
        )

def test_fidelity_fail_inference():
    state_vector = cirq.one_hot(shape=(4,), dtype=np.complex128)
    state_tensor = np.reshape(state_vector, (2, 2))
    with pytest.raises(ValueError, match='Please specify'):
        _ = cirq.fidelity(state_tensor, 4)

def test_fidelity_bad_shape():
    with pytest.raises(ValueError, match='Invalid quantum state'):
        _ = cirq.fidelity(np.array([[[1.0]]]), np.array([[[1.0]]]), qid_shape=(1,))

def test_von_neumann_entropy():
    # 1x1 matrix
    assert cirq.von_neumann_entropy(np.array([[1]])) == 0
    # An EPR pair state (|00> + |11>)(<00| + <11|)
    assert (
        cirq.von_neumann_entropy(0.5 * np.array([1, 0, 0, 1] * np.array([[1], [0], [0], [1]]))) == 0
    )
    # Maximally mixed state
    # yapf: disable
    assert cirq.von_neumann_entropy(np.array(
        [[0.5, 0],
        [0, 0.5]])) == 1
    # yapf: enable
    # 2x2 random unitary, each column as a ket, each ket as a density matrix,
    # linear combination of the two with coefficients 0.1 and 0.9
    res = cirq.testing.random_unitary(2)
    first_column = res[:, 0]
    first_density_matrix = 0.1 * np.outer(first_column, np.conj(first_column))
    second_column = res[:, 1]
    second_density_matrix = 0.9 * np.outer(second_column, np.conj(second_column))
    assert np.isclose(
        cirq.von_neumann_entropy(first_density_matrix + second_density_matrix), 0.4689, atol=1e-04
    )

    assert np.isclose(
        cirq.von_neumann_entropy(np.diag([0, 0, 0.1, 0, 0.2, 0.3, 0.4, 0])), 1.8464, atol=1e-04
    )
    # Random NxN matrix
    probs = np.random.exponential(size=N)
    probs /= np.sum(probs)
    mat = U @ (probs * U).T.conj()

    np.testing.assert_allclose(
        cirq.von_neumann_entropy(mat), -np.sum(probs * np.log(probs) / np.log(2))
    )
    # QuantumState object
    assert (
        cirq.von_neumann_entropy(cirq.quantum_state(np.array([[0.5, 0], [0, 0.5]]), qid_shape=(2,)))
        == 1
    )
    assert (
        cirq.von_neumann_entropy(
            cirq.quantum_state(np.array([[0.5, 0.5], [0.5, 0.5]]), qid_shape=(2, 2))
        )
        == 0
    )

def test_deprecated():
    density_matrix = np.array([[0.5, 0], [0, 0.5]])
    with cirq.testing.assert_deprecated('density_matrix', 'state', deadline="v0.11"):
        # pylint: disable=unexpected-keyword-arg,no-value-for-parameter
        _ = cirq.von_neumann_entropy(density_matrix=density_matrix)

max_stars_repo_name: JiahaoYao/Cirq

max_issues_repo_name: JiahaoYao/Cirq

max_forks_repo_name: JiahaoYao/Cirq


text:
import numpy as np
import torch
from torch import nn
import torch.nn.functional as F
from gym import make

np.random.seed(1)
torch.manual_seed(1)

# Superparameters
MAX_EPISODE = 3000
DISPLAY_REWARD_THRESHOLD = 200  # renders environment if total episode reward is greater then this threshold
MAX_EP_STEPS = 1000  # maximum time step in one episode
RENDER = True  # rendering wastes time
GAMMA = 0.9  # reward discount in TD error
LR_A = 0.001  # learning rate for actor
LR_C = 0.01  # learning rate for critic

env = make('CartPole-v1')
env.seed(1)
env = env.unwrapped

print("env.action_space :", env.action_space)
print("env.observation_space :", env.observation_space)
print("env.observation_space.high :", env.observation_space.high)
print("env.observation_space.low :", env.observation_space.low)

N_F = env.observation_space.shape[0]
N_A = env.action_space.n

class PolicyNet(nn.Module):
    def __init__(self, n_actions, n_features):
        super(PolicyNet, self).__init__()
        hidden_units = 20
        self.fc_layer = nn.Sequential(nn.Linear(n_features, hidden_units),
                                      nn.Sigmoid(),
                                      nn.Linear(hidden_units, n_actions),
                                      nn.Softmax(dim=-1))

    def forward(self, x):
        output = self.fc_layer(x)
        return output

class CriticNet(nn.Module):
    def __init__(self, n_features):
        super(CriticNet, self).__init__()
        hidden_units = 20
        self.fc_layer = nn.Sequential(nn.Linear(n_features, hidden_units),
                                      nn.ReLU(),
                                      nn.Linear(hidden_units, 1))

    def forward(self, x):
        output = self.fc_layer(x)
        return output

class Actor(object):
    def __init__(self, n_features, n_actions, lr=0.001):
        self.actor_net = PolicyNet(n_actions, n_features)
        self.n_features = n_features
        self.n_actions = n_actions
        self.lr = lr

        self.optimizer = torch.optim.Adam(self.actor_net.parameters(),
                                          self.lr)
        self.cost_his = []

    def learn(self, s, a, td):
        state = torch.Tensor(s[np.newaxis, :])
        torch_acts = torch.as_tensor(a)
        torch_acts_one_hot = F.one_hot(torch_acts, num_classes=self.n_actions)
        torch_td_error = torch.Tensor(td).reshape(-1, 1).detach()
        all_act_prob = self.actor_net(state)

        exp_v = torch.log(all_act_prob) * torch_acts_one_hot * torch_td_error
        loss = torch.mean(-exp_v)
        self.optimizer.zero_grad()
        loss.backward()
        self.optimizer.step()
        self.cost_his.append(loss.data.numpy())
        return exp_v

    def choose_action(self, observation):
        state = torch.Tensor(observation[np.newaxis, :])
        prob_weights = self.actor_net(state)
        action_idx = prob_weights.reshape(-1, ).multinomial(num_samples=1).numpy()[0]
        return action_idx

class Critic(object):
    def __init__(self, n_features, lr=0.01):
        self.critic_net = CriticNet(n_features)
        self.n_features = n_features
        self.lr = lr
        self.optimizer = torch.optim.Adam(self.critic_net.parameters(),
                                          self.lr)
        self.cost_his = []
        self.loss_function = torch.nn.MSELoss()

    def learn(self, s, r, s_):
        s, s_ = torch.Tensor(s[np.newaxis, :]), torch.Tensor(s_[np.newaxis, :])

        v = self.critic_net(s)
        v_ = self.critic_net(s_).detach()
        td_error = r + GAMMA * v_ - v
        loss = self.loss_function(v, r + GAMMA * v_)
        # loss = torch.mean(torch.square(td_error))

        self.optimizer.zero_grad()
        loss.backward()
        self.optimizer.step()
        self.cost_his.append(loss.data.numpy())

        return td_error

actor = Actor(n_features=N_F, n_actions=N_A, lr=LR_A)
critic = Critic(n_features=N_F,
                lr=LR_C)  # we need a good teacher, so the teacher should learn faster than the actor

for i_episode in range(MAX_EPISODE):
    s = env.reset()
    t = 0
    track_r = []
    while True:
        if RENDER:
            env.render()

        a = actor.choose_action(s)

        s_, r, done, info = env.step(a)

        if done:
            r = -20

        track_r.append(r)

        td_error = critic.learn(s, r, s_)  # gradient = grad[r + gamma * V(s_) - V(s)]
        actor.learn(s, a, td_error)  # true_gradient = grad[logPi(s,a) * td_error]

        s = s_
        t += 1

        if done or t >= MAX_EP_STEPS:
            ep_rs_sum = sum(track_r)

            if 'running_reward' not in globals():
                running_reward = ep_rs_sum
            else:
                running_reward = running_reward * 0.95 + ep_rs_sum * 0.05
            if running_reward > DISPLAY_REWARD_THRESHOLD: RENDER = True  # rendering
            print("episode:", i_episode, "  reward:", int(running_reward))
            break


text:
"""export classical feature extractor (not tested)
"""

import argparse
import time
import csv
import yaml
import os
import logging
from pathlib import Path
import torch
import cv2

import numpy as np
from tqdm import tqdm

from tensorboardX import SummaryWriter

from utils.utils import tensor2array, save_checkpoint, load_checkpoint, save_path_formatter
from settings import EXPER_PATH
from utils.loader import dataLoader, modelLoader, pretrainedLoader
from utils.utils import getWriterPath

# from utils.logging import *

def export_descriptor(config, output_dir, args):
    '''
    1) input 2 images, output keypoints and correspondence

    :param config:
    :param output_dir:
    :param args:
    :return:
    '''
    # config
    # device = torch.device("cpu")
    device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu")

    logging.info('train on device: %s', device)
    with open(os.path.join(output_dir, 'config.yml'), 'w') as f:
        yaml.dump(config, f, default_flow_style=False)
    writer = SummaryWriter(getWriterPath(task=args.command, date=True))

    ## save data
    from pathlib import Path
    # save_path = save_path_formatter(config, output_dir)
    save_path = Path(output_dir)
    save_output = save_path
    save_output = save_output / 'predictions'
    save_path = save_path / 'checkpoints'
    logging.info('=> will save everything to {}'.format(save_path))
    os.makedirs(save_path, exist_ok=True)
    os.makedirs(save_output, exist_ok=True)

    # data loading
    from utils.loader import dataLoader_test as dataLoader
    data = dataLoader(config, dataset='hpatches')
    test_set, test_loader = data['test_set'], data['test_loader']

    from utils.print_tool import datasize
    datasize(test_loader, config, tag='test')

    from imageio import imread
    def load_as_float(path):
        return imread(path).astype(np.float32) / 255

    def squeezeToNumpy(tensor_arr):
        return tensor_arr.detach().cpu().numpy().squeeze()

    outputMatches = True
    count = 0
    max_length = 5
    method = config['model']['method']
    # tracker = PointTracker(max_length, nn_thresh=fe.nn_thresh)

    # for sample in tqdm(enumerate(test_loader)):
    for i, sample in tqdm(enumerate(test_loader)):

        img_0, img_1 = sample['image'], sample['warped_image']

        imgs_np, imgs_fil = [], []
        # first image, no matches
        imgs_np.append(img_0.numpy().squeeze())
        imgs_np.append(img_1.numpy().squeeze())
        # H, W = img.shape[1], img.shape[2]
        # img = img.view(1,1,H,W)

        ##### add opencv functions here #####
        def classicalDetectors(image, method='sift'):
            """
            # sift keyframe detectors and descriptors
            """
            image = image*255
            round_method = False
            if round_method == True:
                from models.classical_detectors_descriptors import classical_detector_descriptor # with quantization
                points, desc = classical_detector_descriptor(image, **{'method': method})
                y, x = np.where(points)
                # pnts = np.stack((y, x), axis=1)
                pnts = np.stack((x, y), axis=1) # should be (x, y)
                ## collect descriptros
                desc = desc[y, x, :]
            else:
                # sift with subpixel accuracy
                from models.classical_detectors_descriptors import SIFT_det as classical_detector_descriptor
                pnts, desc = classical_detector_descriptor(image, image)

            print("desc shape: ", desc.shape)
            return pnts, desc

        pts_list = []
        pts, desc_1 = classicalDetectors(imgs_np[0], method=method)
        pts_list.append(pts)
        print("total points: ", pts.shape)
        '''
        pts: list [numpy (N, 2)]
        desc: list [numpy (N, 128)]
        '''
        # save keypoints
        pred = {}
        pred.update({
            # 'image_fil': imgs_fil[0],
            'image': imgs_np[0],
        })
        pred.update({'prob': pts,
                     'desc': desc_1})

        # second image, output matches

        pred.update({
            'warped_image': imgs_np[1],
            # 'warped_image_fil': imgs_fil[1],
        })
        pts, desc_2 = classicalDetectors(imgs_np[1], method=method)
        pts_list.append(pts)

        # if outputMatches == True:
        #     tracker.update(pts, desc)
        # pred.update({'matches': matches.transpose()})

        print("total points: ", pts.shape)
        pred.update({'warped_prob': pts,
                     'warped_desc': desc_2,
                     'homography': squeezeToNumpy(sample['homography'])
                     })

        ## get matches
        data = get_sift_match(sift_kps_ii=pts_list[0], sift_des_ii=desc_1, sift_kps_jj=pts_list[1], sift_des_jj=desc_2, if_BF_matcher=True)
        matches = data['match_quality_good']
        print(f"matches: {matches.shape}")
        matches_all = data['match_quality_all']
        pred.update({
            'matches': matches,
            'matches_all': matches_all
        })
        # clean last descriptor
        '''
        pred:
            'image': np(320,240)
            'prob' (keypoints): np (N1, 2)
            'desc': np (N2, 256)
            'warped_image': np(320,240)
            'warped_prob' (keypoints): np (N2, 2)
            'warped_desc': np (N2, 256)
            'homography': np (3,3)
            'matches': np (N3, 4)
        '''

        # save data
        from pathlib import Path
        filename = str(count)
        path = Path(save_output, '{}.npz'.format(filename))
        np.savez_compressed(path, **pred)
        count += 1

    print("output pairs: ", count)
    save_file = save_output / "export.txt"
    with open(save_file, "a") as myfile:
        myfile.write("output pairs: " + str(count) + '\n')
    pass

def get_sift_match(sift_kps_ii, sift_des_ii, sift_kps_jj, sift_des_jj, if_BF_matcher=True):
    # select which kind of matcher
    if (
        if_BF_matcher
    ):  # OpenCV sift matcher must be created inside each thread (because it does not support sharing across threads!)
        bf = cv2.BFMatcher(normType=cv2.NORM_L2)
        sift_matcher = bf
    else:
        FLANN_INDEX_KDTREE = 0
        index_params = dict(algorithm=FLANN_INDEX_KDTREE, trees=5)
        search_params = dict(checks=50)
        flann = cv2.FlannBasedMatcher(index_params, search_params)
        sift_matcher = flann

    all_ij, good_ij, quality_good, quality_all, cv_matches = get_sift_match_idx_pair(
        sift_matcher, sift_des_ii.copy(), sift_des_jj.copy()
    )
    if all_ij is None:
        logging.warning(
            "KNN match failed dumping %s frame %d-%d. Skipping" % (dump_dir, ii, jj)
        )
    # dump_ij_idx_file = dump_dir / "ij_idx_{}-{}".format(ii, jj)
    # dump_ij_quality_file = dump_dir / "ij_quality_{}-{}".format(ii, jj)
    # dump_ij_match_quality_file = dump_dir / "ij_match_quality_{}-{}".format(ii, jj)

    # print(good_ij, good_ij.shape)
    match_quality_good = np.hstack(
        (sift_kps_ii[good_ij[:, 0]], sift_kps_jj[good_ij[:, 1]], quality_good)
    )  # [[x1, y1, x2, y2, dist_good, ratio_good]]
    match_quality_all = np.hstack(
        (sift_kps_ii[all_ij[:, 0]], sift_kps_jj[all_ij[:, 1]], quality_all)
    )  # [[x1, y1, x2, y2, dist_good, ratio_good]]
    return {'match_quality_good': match_quality_good,
    'match_quality_all': match_quality_all, 'cv_matches': cv_matches}

def get_sift_match_idx_pair(sift_matcher, des1, des2):
    """
    do matchings, test the quality of matchings
    """
    try:
        matches = sift_matcher.knnMatch(
            des1, des2, k=2
        )  # another option is https://github.com/MagicLeapResearch/SuperPointPretrainedNetwork/blob/master/demo_superpoint.py#L309
    except Exception as e:
        logging.error(traceback.format_exception(*sys.exc_info()))
        return None, None
    # store all the good matches as per Lowe's ratio test.
    good = []
    all_m = []
    quality_good = []
    quality_all = []
    for m, n in matches:
        all_m.append(m)
        if m.distance < 0.8 * n.distance:
            good.append(m)
            quality_good.append([m.distance, m.distance / n.distance])
        quality_all.append([m.distance, m.distance / n.distance])

    good_ij = [[mat.queryIdx for mat in good], [mat.trainIdx for mat in good]]
    all_ij = [[mat.queryIdx for mat in all_m], [mat.trainIdx for mat in all_m]]
    return (
        np.asarray(all_ij, dtype=np.int32).T.copy(),
        np.asarray(good_ij, dtype=np.int32).T.copy(),
        np.asarray(quality_good, dtype=np.float32).copy(),
        np.asarray(quality_all, dtype=np.float32).copy(),
        matches
    )

if __name__ == '__main__':
    # global var
    torch.set_default_tensor_type(torch.FloatTensor)
    logging.basicConfig(format='[%(asctime)s %(levelname)s] %(message)s',
                        datefmt='%m/%d/%Y %H:%M:%S', level=logging.INFO)

    # add parser
    parser = argparse.ArgumentParser()
    subparsers = parser.add_subparsers(dest='command')

    # export command
    p_train = subparsers.add_parser('export_descriptor')
    p_train.add_argument('config', type=str)
    p_train.add_argument('exper_name', type=str)
    # p_train.add_argument('exper', type=str)
    p_train.add_argument('--correspondence', action='store_true')
    p_train.add_argument('--eval', action='store_true')
    p_train.add_argument('--debug', action='store_true', default=False,
                         help='turn on debuging mode')
    p_train.set_defaults(func=export_descriptor)

    args = parser.parse_args()
    with open(args.config, 'r') as f:
        config = yaml.load(f)

    output_dir = os.path.join(EXPER_PATH, args.exper_name)
    os.makedirs(output_dir, exist_ok=True)

    # with capture_outputs(os.path.join(output_dir, 'log')):
    logging.info('Running command {}'.format(args.command.upper()))
    args.func(config, output_dir, args)

    # global variables

    # main()

max_stars_repo_name: Dai-z/pytorch-superpoint

max_issues_repo_name: Dai-z/pytorch-superpoint

max_forks_repo_name: Dai-z/pytorch-superpoint


text:
import numpy as np

from tools import getLink

class linkDensity:
    #dt is the threshold to decide whether the domTree(sub-tree) belongs to one cluster
    dt = 0
    #init the linkDensity
    def __init__(self, domTree, alpha = 1, beta = 1):
        self.alpha = alpha
        self.beta = beta
        self.dt = self.denThres(domTree)

    #get the density threshold for judging whether this dom-node is one cluster
    #bigger value means higher possibility to be one cluster
    #   input:  dom-tree
    #   output: the density threshold
    def denThres(self, domTree):
        linkList = getLink(domTree)
        linkStrList = []
        for link in linkList:
            linkStrList.append(link.text)
            imgList = link.find_all('img')
            if not imgList:
                linkStrList.append('string')
            for img in link.find_all('img'):
                if 'alt' in img.attrs:
                    linkStrList.append(img['alt'])
                else:
                    linkStrList.append('string')
        linkStrLen = sum([len(linkStr.split()) for linkStr in linkStrList])
        dt = self.alpha*(linkStrLen+self.beta)/(float(len(domTree.text.split()))+linkStrLen+self.beta)
        return dt

    #get the link density in a dom-tree
    #   input: one dom-tree
    #   output: the link density value, if any error return None
    def getLinkDen(self, subDomList):
        linkList = []
        for subDom in subDomList:
            linkList += getLink(subDom)
        linkStrList = []
        for link in linkList:
            linkStrList.append(link.text)
            imgList = link.find_all('img')
            if not imgList:
                linkStrList.append('string')
            for img in link.find_all('img'):
                if 'alt' in img.attrs:
                    linkStrList.append(img['alt'])
                else:
                    linkStrList.append('string')
        linkStrLen = sum([len(linkStr.split()) for linkStr in linkStrList])
        ld = (linkStrLen+self.beta)/(float(len(subDom.text.split()))+linkStrLen+self.beta)
        return ld

    #decide whether links in a sub Dom-tree are in one block
    def isOneCluster(self, subDom):
        if self.alpha == 0:
            return True
        return self.getLinkDen(subDom)>self.dt

max_stars_repo_name: zhaokv/Navigation-Objects-Extraction

max_issues_repo_name: zhaokv/Navigation-Objects-Extraction

max_forks_repo_name: zhaokv/Navigation-Objects-Extraction


text:
# set up environment
# define command line parameters
# define location of input data
# create output directories
# start the class FeatureCorrelations

import argparse
import os

#import matplotlib
#matplotlib.use("Agg")

import matplotlib.pyplot as plt
import matplotlib.cm as cm
import csv
import pathlib
import numpy as np
import pandas as pd
import seaborn as sns
from datetime import datetime
from sklearn.cluster import DBSCAN
from sklearn import metrics
from sklearn.preprocessing import StandardScaler

###############################################################################
#
#  define command line arguments
#
###############################################################################

def parse_command_line():
  '''defining the command line input to make it runable'''
  parser = argparse.ArgumentParser(description='various plots for feature analysis')

  parser.add_argument(
    '--input',
    type=str,
    dest="input",
    default="",
    help='The input CSV file')

  parser.add_argument(
    '--outdir',
    type=str,
    dest='outdir',
    default='',
    help='Specify output directory')

  args = parser.parse_args()
  if args.input == '':
    parser.print_help()
    exit(0)
  return args

###############################################################################
#
#  load the data from CSV file and creating output directory
#
###############################################################################

def load_data(csv_path):
  '''load the raw data as stored in CSV file'''
  return pd.read_csv(csv_path, na_filter=False, skipinitialspace=True, thousands=',')

def make_output_folder(outdir):
  name = os.path.join(outdir, 'dbscan_clustering')
  os.makedirs(name, exist_ok=True)
  return name

###############################################################################
#
#  class to analyse correlations between features
#
###############################################################################

class DBscan(object):
  '''A class to help analyse the data;
  try to identify linear correlations in the data;
  calculate Pearson Correlation Coefficient with and without
  p-values; create a scatter matrix; inout data must not contain
  any strings or NaN values; also remove any columns with
  categorical data or transform them first; remove any text except column labels'''

  def __init__(self, data, dbscan_clustering):
    self.data = data
    self.dbscan_clustering = dbscan_clustering
    self.prepare_metrix_data()
    self.clustering()

  ###############################################################################
  #
  #  creating 3 data frames specific to the three development milestones I had
  #  1--> directly from data processing
  #  2--> after adding protein information
  #  3--> carrying out some further column transformations
  #
  ###############################################################################

  def prepare_metrix_data(self):
    '''Function to create smaller dataframes for directly after dataprocessing, after
       adding some protein information and after carrying out some custom solumn
       transformations.
    ******
    Input: large data frame
    Output: smaller dataframes; database, man_add, transform
    '''
    print('*' *80)
    print('*    Preparing input dataframe')
    print('*' *80)

    self.X_data = self.data[["mapCC", "res_frag_ratio"]]

    self.X_data = self.X_data.fillna(0)

    self.X_data_scaled = StandardScaler().fit_transform(self.X_data)

    with open(os.path.join(self.dbscan_clustering, 'dbscan_clustering.txt'), 'a') as text_file:
      text_file.write('Created the following dataframe: X_data \n')
      #text_file.write(str(self.X_data.columns)+'\n')

    #print(self.X_data_scaled.columns)

################################################################################
#
#  plotting probability density function for each feature to evaluate continous
#  data and its cumulative distribution function
#
################################################################################

  def clustering(self):
    '''Cluster samples to find classes.'''
    print('*' *80)
    print('*    Clustering samples')
    print('*' *80)

    datestring = datetime.strftime(datetime.now(), '%Y%m%d_%H%M')

    db = DBSCAN(eps=0.4, min_samples=5).fit(self.X_data_scaled)
    core_samples_mask = np.zeros_like(db.labels_, dtype=bool)
    core_samples_mask[db.core_sample_indices_] = True
    labels = db.labels_

    # Number of clusters in labels, ignoring noise if present.
    n_clusters_ = len(set(labels)) - (1 if -1 in labels else 0)
    n_noise_ = list(labels).count(-1)

    print('Estimated number of clusters: %d' % n_clusters_)
    print('Estimated number of noise points: %d' % n_noise_)
#    print("Homogeneity: %0.3f" % metrics.homogeneity_score(labels_true, labels))
#    print("Completeness: %0.3f" % metrics.completeness_score(labels_true, labels))
#    print("V-measure: %0.3f" % metrics.v_measure_score(labels_true, labels))
#    print("Adjusted Rand Index: %0.3f"
#          % metrics.adjusted_rand_score(labels_true, labels))
#    print("Adjusted Mutual Information: %0.3f"
#          % metrics.adjusted_mutual_info_score(labels_true, labels))
#    print("Silhouette Coefficient: %0.3f"
#          % metrics.silhouette_score(self.X_data_scaled, labels))

    # #############################################################################
    # Plot result
    import matplotlib.pyplot as plt

    # Black removed and is used for noise instead.
    unique_labels = set(labels)
    colors = [plt.cm.Spectral(each)
              for each in np.linspace(0, 1, len(unique_labels))]
    for k, col in zip(unique_labels, colors):
        if k == -1:
            # Black used for noise.
            col = [0, 0, 0, 1]

        class_member_mask = (labels == k)

        xy = self.X_data_scaled[class_member_mask & core_samples_mask]
        plt.plot(xy[:, 0], xy[:, 1], 'o', markerfacecolor=tuple(col),
                 markeredgecolor='k', markersize=14)

        xy = self.X_data_scaled[class_member_mask & ~core_samples_mask]
        plt.plot(xy[:, 0], xy[:, 1], 'o', markerfacecolor=tuple(col),
                 markeredgecolor='k', markersize=6)

    plt.title('Estimated number of clusters: %d' % n_clusters_)
    plt.savefig(os.path.join(self.dbscan_clustering, 'alldata_clusters_'+datestring+'.png'))
    plt.show()

#        plt.savefig(os.path.join(self.k_means_clustering, 'Silhouette_autosharp_clusters_'+str(n_clusters)+'_'+datestring+'.png'))
#        plt.close()
    #plt.show()

def run():
  args = parse_command_line()

  ###############################################################################

  #look at the imported data to get an idea what we are working with
  data = load_data(args.input)

  dbscan_clustering = make_output_folder(args.outdir)

  ###############################################################################

  dbscan_clustering = DBscan(data, dbscan_clustering)


text:
# coding: utf-8
import io
import json
import os
from pathlib import PurePath

import fire
import matplotlib.pyplot as plt
import numpy as np
from seaborn import heatmap

def load(filename):
    replacement_matrix = {}
    with open(filename) as f:
        key_names = f.readline().strip().split()
        for line in f:
            key_values = line.strip().split()
            key, values = key_values[0], key_values[1:]
            assert len(key_names) == len(values), "\nkeys(%d): %s\nvalues(%d):%s" % (
                len(key_names), key_names, len(values), values
            )
            for cmp_key, value in zip(key_names, values):
                if key not in replacement_matrix:
                    replacement_matrix[key] = {}
                replacement_matrix[key][cmp_key] = int(value)
    return replacement_matrix

class BLOSUM(object):
    def __init__(self, matrix_type='62'):
        self.matrix = load(os.path.join(PurePath(__file__).parents[1], ("meta/BLOSUM%s" % matrix_type)))

    def __call__(self, key, cmp_key):
        return self.matrix[key][cmp_key]

    def __str__(self):
        matrix = [[""] * len(self.matrix) for _ in range(len(self.matrix))]
        sorted_keys = sorted(self.matrix.keys())
        for i, key in enumerate(sorted_keys):
            for j, cmp_key in enumerate(sorted_keys):
                matrix[i][j] = "%2s" % str(self(key, cmp_key))

        ret = []
        ret.append("\t".join([" "] + ["%2s" % key for key in sorted_keys]))
        for i, key in enumerate(sorted_keys):
            ret.append("\t".join([key] + matrix[i]))
        return "\n".join(ret)

    def __repr__(self):
        return str(self)

class SequenceComparePlayer(object):
    def __init__(self, matrix_type="62"):
        self.blosum = BLOSUM(matrix_type)

    def __call__(self, seq1: str, seq2: str) -> int:
        assert len(seq1) == len(seq2)
        score = 0
        if seq1 == seq2:
            seq1 = seq1.replace("-", "")
            seq2 = seq2.replace("-", "")
        for key1, key2 in zip(seq1, seq2):
            if key1 == "-" or key2 == "-":
                score += 11
            else:
                score += self.blosum(key1, key2)
        return score

def run(src, tar=None, itar=None, viz=True):
    seqs = []
    with open(src) as f:
        for line in f:
            line = line.strip()
            if line:
                seqs.append(line)

    result = [[0] * len(seqs) for _ in range(len(seqs))]

    scp = SequenceComparePlayer()

    for i, seq1 in enumerate(seqs):
        for j, seq2 in enumerate(seqs):
            result[i][j] = scp(seq1, seq2)

    if viz or itar:
        viz_result = np.array(result)
        plt.clf()
        heatmap(viz_result)
        if itar:
            plt.savefig(itar)
        if viz:
            plt.show()

    if tar:
        with io.open(tar, "w") as wf:
            json.dump(result, wf)

if __name__ == '__main__':
    fire.Fire(run)

max_stars_repo_name: tswsxk/seqcp

max_issues_repo_name: tswsxk/seqcp

max_forks_repo_name: tswsxk/seqcp


text:
from numpy import array

from automaton import Neighborhood
from automaton.CellState import CellState, cell_state_from

def get_processing_function_bundle(dimensions, neighborhood_type, birth_nums, survival_nums):
    neighborhood = None
    processing_type = dimensions + '' + neighborhood_type
    if processing_type == '2A':
        neighborhood = Neighborhood.moore2d
    elif processing_type == '2B':
        neighborhood = Neighborhood.neumann2d

    def processing_function(world, x, y, extended_output=False):
        output = neighborhood(world, x, y, extended_output)
        n = output[- 1]
        if world.get_cell(x, y) == CellState.Alive:
            if n in survival_nums:
                new_state = CellState.Alive
            else:
                new_state = CellState.Dead
        else:
            if n in birth_nums:
                new_state = CellState.Alive
            else:
                new_state = CellState.Dead
        if extended_output:
            output[len(output) - 1] = new_state
            return output
        else:
            return new_state

    return processing_function, neighborhood

def get_neural_processing_function_bundle(neighborhood, neural_network):
    def processing_function(world, x, y, extended_output=False):
        output = neighborhood(world, x, y, True)
        ret = neural_network.think_output(array([float(i) for i in output[:-1]]))
        if ret > 0.5:
            ret = 1
        else:
            ret = 0
        new_state = cell_state_from(ret)
        if extended_output:
            output[len(output) - 1] = new_state
            return output
        else:
            return new_state

    return processing_function, neighborhood

reason: from numpy


text:
'''
#                           This is an 80 character line                       #
Purpose: read data from multiple text files and compare (plot) in one convenient
image :) Wow, go you
'''

import sys
import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt

# Use this to back-calculate time in units of Brownian tau
kT = 1.0                        # temperature
threeEtaPiSigma = 1.0           # drag coefficient
sigma = 1.0                     # particle diameter
D_t = kT / threeEtaPiSigma      # translational diffusion constant
D_r = (3.0 * D_t) / (sigma**2)  # rotational diffusion constant
tauBrown = (sigma**2) / D_t     # brownian time scale (invariant)

# Define what functions you'll need here
def getFromTxt(fname, first, last):
    "Takes a string, text before and after desired text, outs text between"
    start = fname.index( first ) + len( first )
    end = fname.index( last, start )
    myTxt = fname[start:end]
    return float(myTxt)
# Above function kindly provided by user "cji" on stackoverflow
#https://stackoverflow.com/questions/3368969/find-string-between-two-substrings

def computeVel(activity):
    "Given particle activity, output intrinsic swim speed"
    velocity = (activity * sigma) / (3 * (1/D_r))
    return velocity

def computeActiveForce(velocity):
    "Given particle activity, output repulsion well depth"
    activeForce = velocity * threeEtaPiSigma
    return activeForce

def computeEps(activeForce):
    "Given particle activity, output repulsion well depth"
    epsilon = activeForce * sigma / 24.0
    return epsilon

def computeTauLJ(epsilon):
    "Given epsilon, compute lennard-jones time unit"
    tauLJ = ((sigma**2) * threeEtaPiSigma) / epsilon
    return tauLJ

# Grab the textfiles to run on
txtFiles = sys.argv[2:] # pass starting at 2 to avoid script path

# Parse out the activities and fractions from the filenames
peA = np.zeros_like(txtFiles, dtype=np.float64)
peB = np.zeros_like(txtFiles, dtype=np.float64)
xA = np.zeros_like(txtFiles, dtype=np.float64)
ep = np.zeros_like(txtFiles, dtype=np.float64)

# This grabs the parameters of each text file
for i in range(0, len(txtFiles)):
    peA[i] = getFromTxt(txtFiles[i], "pa", "_pb")
    peB[i] = getFromTxt(txtFiles[i], "pb", "_xa")
    xA[i] = getFromTxt(txtFiles[i], "xa", "_ep")
    ep[i] = getFromTxt(txtFiles[i], "ep", ".txt")

partFracA = xA/100.0

# Initialize the arrays so that all text files fit in each one
tsteps = np.zeros(len(txtFiles), dtype=np.ndarray)
gasA = np.zeros(len(txtFiles), dtype=np.ndarray)
gasB = np.zeros(len(txtFiles), dtype=np.ndarray)
gasAll = np.zeros(len(txtFiles), dtype=np.ndarray)
denseA = np.zeros(len(txtFiles), dtype=np.ndarray)
denseB = np.zeros(len(txtFiles), dtype=np.ndarray)
denseAll = np.zeros(len(txtFiles), dtype=np.ndarray)
lgClust = np.zeros(len(txtFiles), dtype=np.ndarray)
MCS = np.zeros(len(txtFiles), dtype=np.ndarray)
sigALL = np.zeros(len(txtFiles), dtype=np.ndarray)
sigAA = np.zeros(len(txtFiles), dtype=np.ndarray)
sigAB = np.zeros(len(txtFiles), dtype=np.ndarray)
sigBB = np.zeros(len(txtFiles), dtype=np.ndarray)
phiEff = np.zeros(len(txtFiles), dtype=np.ndarray)
lgClustA = np.zeros(len(txtFiles), dtype=np.ndarray)
totClustA = np.zeros(len(txtFiles), dtype=np.ndarray)
lgClustD = np.zeros(len(txtFiles), dtype=np.ndarray)
dpDensity = np.zeros(len(txtFiles), dtype=np.ndarray)
gpDensity = np.zeros(len(txtFiles), dtype=np.ndarray)
dt = np.zeros(len(txtFiles), dtype=np.float64)
# It is absolutley archaic that I have to initialize these like this ^

# Pull data from text files
for i in range(0, len(txtFiles)):
    tsteps[i],\
    gasA[i],\
    gasB[i],\
    gasAll[i],\
    denseA[i],\
    denseB[i],\
    denseAll[i],\
    lgClust[i],\
    MCS[i],\
    sigALL[i],\
    sigAA[i],\
    sigAB[i],\
    sigBB[i],\
    phiEff[i],\
    lgClustA[i],\
    totClustA[i],\
    lgClustD[i],\
    dpDensity[i],\
    gpDensity[i] = np.loadtxt(txtFiles[i], skiprows=1, unpack=True)

    partNum = gasAll[i][0]
    gasA[i] /= partNum * partFracA[i]
    gasB[i] /= partNum * (1-partFracA[i])
    gasAll[i] /= partNum
    denseA[i] /= partNum * partFracA[i]
    denseB[i] /= partNum * (1-partFracA[i])
    denseAll[i] /= partNum
    lgClust[i] /= partNum

# Compute Brownian time
for i in range(0, len(txtFiles)):
    # Compute parameters from activities
    if peA[i] != 0:                     # A particles are NOT Brownian
        vA = computeVel(peA[i])
        FpA = computeActiveForce(vA)
        epsA = computeEps(FpA)
        tauA = computeTauLJ(epsA)
    else:                               # A particles are Brownian
        vA = 0.0
        FpA = 0.0
        epsA = kT
        tauA = computeTauLJ(epsA)

    if peB[i] != 0:                     # B particles are NOT Brownian
        vB = computeVel(peB[i])
        FpB = computeActiveForce(vB)
        epsB = computeEps(FpB)
        tauB = computeTauLJ(epsB)
    else:                               # B particles are Brownian
        vB = 0.0
        FpB = 0.0
        epsB = kT
        tauB = computeTauLJ(epsB)

    # Get adjusted time units
    tauLJ = (tauA if (tauA <= tauB) else tauB)  # use the smaller tauLJ
    ratio = tauLJ / tauBrown                    # get rid of LJ units
    dt[i] = ratio * 0.00001                     # tstep size
    tsteps[i] *= dt[i]

# Reorder things in order of increasing timestep
for i in range(0, len(txtFiles)):
    for j in range(0, len(tsteps[i])):
        for k in range(0, len(tsteps[i])):
            # Need to swap
            if tsteps[i][j] < tsteps[i][k] and j > k:
                tmp = tsteps[i][j]
                tsteps[i][j] = tsteps[i][k]
                tsteps[i][k] = tmp
                tmp = gasA[i][j]
                gasA[i][j] = gasA[i][k]
                gasA[i][k] = tmp
                tmp = gasB[i][j]
                gasB[i][j] = gasB[i][k]
                gasB[i][k] = tmp
                tmp = gasAll[i][j]
                gasAll[i][j] = gasAll[i][k]
                gasAll[i][k] = tmp
                tmp = denseA[i][j]
                denseA[i][j] = denseA[i][k]
                denseA[i][k] = tmp
                tmp = denseB[i][j]
                denseB[i][j] = denseB[i][k]
                denseB[i][k] = tmp
                tmp = denseAll[i][j]
                denseAll[i][j] = denseAll[i][k]
                denseAll[i][k] = tmp
                tmp = lgClust[i][j]
                lgClust[i][j] = lgClust[i][k]
                lgClust[i][k] = tmp
                tmp = MCS[i][j]
                MCS[i][j] = MCS[i][k]
                MCS[i][k] = tmp
                tmp = sigALL[i][j]
                sigALL[i][j] = sigALL[i][k]
                sigALL[i][k] = tmp
                tmp = sigAA[i][j]
                sigAA[i][j] = sigAA[i][k]
                sigAA[i][k] = tmp
                tmp = sigAB[i][j]
                sigAB[i][j] = sigAB[i][k]
                sigAB[i][k] = tmp
                tmp = sigBB[i][j]
                sigBB[i][j] = sigBB[i][k]
                sigBB[i][k] = tmp
                tmp = phiEff[i][j]
                phiEff[i][j] = phiEff[i][k]
                phiEff[i][k] = tmp
                tmp = lgClustA[i][j]
                lgClustA[i][j] = lgClustA[i][k]
                lgClustA[i][k] = tmp
                tmp = totClustA[i][j]
                totClustA[i][j] = totClustA[i][k]
                totClustA[i][k] = tmp
                tmp = lgClustD[i][j]
                lgClustD[i][j] = lgClustD[i][k]
                lgClustD[i][k] = tmp
                tmp = dpDensity[i][j]
                dpDensity[i][j] = dpDensity[i][k]
                dpDensity[i][k] = tmp
                tmp = gpDensity[i][j]
                gpDensity[i][j] = gpDensity[i][k]
                gpDensity[i][k] = tmp

# Now plot everything

# Gas phase
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], gasA[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Gas Phase')
plt.ylim(0,1)
plt.savefig('gasPhaseA.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], gasB[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Gas Phase')
plt.ylim(0,1)
plt.savefig('gasPhaseB.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], gasAll[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Gas Phase')
plt.ylim(0,1)
plt.savefig('gasPhaseAll.png', dpi=1000)
plt.close()

# Dense phase
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], denseA[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Dense Phase')
plt.ylim(0,1)
plt.savefig('densePhaseA.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], denseB[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Dense Phase')
plt.ylim(0,1)
plt.savefig('densePhaseB.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], denseAll[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Dense Phase')
plt.ylim(0,1)
plt.savefig('densePhaseAll.png', dpi=1000)
plt.close()

# Largest Cluster
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], lgClust[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Fraction in Largest Cluster')
plt.ylim(0,1)
plt.savefig('largestCluster.png', dpi=1000)
plt.close()

# MCS
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], MCS[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Mean Cluster Size')
plt.savefig('MCS.png', dpi=1000)
plt.close()

# Largest Cluster Area
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], lgClustA[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Largest Cluster Area')
plt.savefig('largest_cluster_area.png', dpi=1000)
plt.close()

# Total Cluster Area
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], totClustA[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Total Cluster Area')
plt.savefig('total_cluster_area.png', dpi=1000)
plt.close()

# Largest Cluster Density
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], lgClustD[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Largest Cluster Density')
plt.savefig('largest_cluster_density.png', dpi=1000)
plt.close()

# Dense phase density
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], dpDensity[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Dense Phase Density')
plt.savefig('dpDensity.png', dpi=1000)
plt.close()

# Gas phase density
for i in range(0, len(txtFiles)):
    plt.plot(tsteps[i], gpDensity[i], label=str(peA[i]))
plt.legend()
plt.xlabel(r'Brownian Time $(\tau_{Brownian})$')
plt.ylabel('Gas Phase Density')
plt.savefig('gpDensity.png', dpi=1000)
plt.close()

################################################################################
# Steady state values
################################################################################
ssGasA = np.zeros((len(txtFiles)), dtype=np.float64)
ssGasB = np.zeros((len(txtFiles)), dtype=np.float64)
ssGasAll = np.zeros((len(txtFiles)), dtype=np.float64)
ssDenseA = np.zeros((len(txtFiles)), dtype=np.float64)
ssDenseB = np.zeros((len(txtFiles)), dtype=np.float64)
ssDenseAll = np.zeros((len(txtFiles)), dtype=np.float64)
ssLgClust = np.zeros((len(txtFiles)), dtype=np.float64)
ssMCS = np.zeros((len(txtFiles)), dtype=np.float64)
ssALL = np.zeros((len(txtFiles)), dtype=np.float64)
ssAA = np.zeros((len(txtFiles)), dtype=np.float64)
ssAB = np.zeros((len(txtFiles)), dtype=np.float64)
ssBB = np.zeros((len(txtFiles)), dtype=np.float64)
ssPhi = np.zeros((len(txtFiles)), dtype=np.float64)
sslgCA = np.zeros((len(txtFiles)), dtype=np.float64)
sstotCA = np.zeros((len(txtFiles)), dtype=np.float64)
sslgCD = np.zeros((len(txtFiles)), dtype=np.float64)
ssDPD = np.zeros((len(txtFiles)), dtype=np.float64)
ssGPD = np.zeros((len(txtFiles)), dtype=np.float64)

# Get steady state values (last 100 dumps)
for i in range(0, len(txtFiles)):
    ssGasA[i] = np.mean(gasA[i][-100:-1])
    ssGasB[i] = np.mean(gasB[i][-100:-1])
    ssGasAll[i] = np.mean(gasAll[i][-100:-1])
    ssDenseA[i] = np.mean(denseA[i][-100:-1])
    ssDenseB[i] = np.mean(denseB[i][-100:-1])
    ssDenseAll[i] = np.mean(denseAll[i][-100:-1])
    ssLgClust[i] = np.mean(lgClust[i][-100:-1])
    ssMCS[i] = np.mean(MCS[i][-100:-1])
    ssALL[i] = np.mean(sigALL[i][-100:-1])
    ssAA[i] = np.mean(sigAA[i][-100:-1])
    ssAB[i] = np.mean(sigAB[i][-100:-1])
    ssBB[i] = np.mean(sigBB[i][-100:-1])
    ssPhi[i] = np.mean(phiEff[i][-100:-1])
    sslgCA[i] = np.mean(lgClustA[i][-100:-1])
    sstotCA[i] = np.mean(totClustA[i][-100:-1])
    sslgCD[i]= np.mean(lgClustD[i][-100:-1])
    ssDPD[i] = np.mean(dpDensity[i][-100:-1])
    ssGPD[i] = np.mean(gpDensity[i][-100:-1])

# Gas Phase #####################

# What is varying? (xA or PeA)
#if xA[0] == xA[1] == 100:
#    plt.scatter(peA, ssGasA, c=peA)
#    plt.xlabel(r'Activity')
#elif xA[0] == xA[1] == 0:
#    plt.scatter(peB, ssGasA, c=peB)
#    plt.xlabel(r'Activity')
#else:
#    ratio = np.zeros_like(peA)
#    ratio = peA / peB
#    plt.scatter(ratio, ssGasA, label=str(ratio), c=ratio)
#    plt.xlabel(r'Activity Ratio')
#
#plt.ylabel('Steady-State Fraction of A-Particles in Gas')
#plt.legend()
##plt.ylim(0,1)
#plt.savefig('SteadyState_gasA.png', dpi=1000)
#plt.close()

for i in range(0, len(txtFiles)):
    # All monodisperse B, use activity
    if xA[i] == xA[i-1] == 0:
        plt.scatter(peB[i], ssGasA[i], c='k')
        plt.xlabel(r'Activity')
    # All monodisperse A, use activity
    elif xA[i] == xA[i-1] == 100:
        plt.scatter(peA[i], ssGasA[i], c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssGasA[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, varying peA
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssGasA[i], c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Fraction of A-Particles in Gas')
plt.savefig('SteadyState_gasA.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssGasB[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssGasB[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssGasB[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssGasB[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Fraction of B-Particles in Gas')
plt.savefig('SteadyState_gasB.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssGasAll[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssGasAll[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssGasAll[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssGasAll[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Fraction of All Particles in Gas')
plt.savefig('SteadyState_gasAll.png', dpi=1000)
plt.close()

# Dense Phase ###################
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssDenseA[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssDenseA[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssDenseA[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssDenseA[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Fraction of A-Particles in Dense Phase')
plt.savefig('SteadyState_denseA.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssDenseB[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssDenseB[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssDenseB[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssDenseB[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Fraction of B-Particles in Dense Phase')
plt.savefig('SteadyState_denseB.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssDenseAll[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssDenseAll[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssDenseAll[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssDenseAll[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Fraction of All Particles in Dense Phase')
plt.savefig('SteadyState_denseAll.png', dpi=1000)
plt.close()

# Largest Cluster ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssLgClust[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssLgClust[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssLgClust[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssLgClust[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Largest Cluster Size')
plt.savefig('SteadyState_lgCluster.png', dpi=1000)
plt.close()

# MCS ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssMCS[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssMCS[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssMCS[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssMCS[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State MCS')
plt.savefig('SteadyState_MCS.png', dpi=1000)
plt.close()

# Diameter All ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        ALL = plt.scatter(peB[i], ssALL[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        ALL = plt.scatter(peA[i], ssALL[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        ALL = plt.scatter(xA[i], ssALL[i], c='k', label='All')
        AA = plt.scatter(xA[i], ssAA[i], c='b', label='AA')
        AB = plt.scatter(xA[i], ssAB[i], c='g', label='AB')
        BB = plt.scatter(xA[i], ssBB[i], c='r', label='BB')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
        legend = plt.legend(bbox_to_anchor=(1.0, 1.0),
                            loc='best',
                            ncol=1,
                            handles=[ALL, AA, AB, BB],
                            labels=['ALL', 'AA', 'AB', 'BB'],
                            framealpha=1.0)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        ALL = plt.scatter(ratio, ssALL[i], c='k', label='All')
        AA = plt.scatter(ratio, ssAA[i], c='b', label='AA')
        AB = plt.scatter(ratio, ssAB[i], c='g', label='AB')
        BB = plt.scatter(ratio, ssBB[i], c='r', label='BB')
        plt.xlabel(r'Activity Ratio')
        legend = plt.legend(bbox_to_anchor=(1.0, 1.0),
                            loc='best',
                            ncol=1,
                            handles=[ALL, AA, AB, BB],
                            labels=['ALL', 'AA', 'AB', 'BB'],
                            framealpha=1.0)

plt.ylabel('Steady-State Diameter')
plt.savefig('SteadyState_diameter.png', dpi=1000)
plt.close()

# Phi Effective ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssPhi[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssPhi[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssPhi[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssPhi[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Phi')
plt.savefig('SteadyState_Phi.png', dpi=1000)
plt.close()

# Steady-state largest cluster area ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], sslgCA[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], sslgCA[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], sslgCA[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, sslgCA[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State largest Cluster area')
plt.savefig('SteadyState_lgCA.png', dpi=1000)
plt.close()

# Steady-state total cluster area ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], sstotCA[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], sstotCA[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], sstotCA[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, sstotCA[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State total Cluster Area')
plt.savefig('SteadyState_totCA.png', dpi=1000)
plt.close()

# Steady-state largest cluster density ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], sslgCD[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], sslgCD[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], sslgCD[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, sslgCD[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State largest cluster density')
plt.savefig('SteadyState_lgCD.png', dpi=1000)
plt.close()

# Cluster Density ###############
for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssDPD[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssDPD[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssDPD[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssDPD[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Dense Phase Density')
plt.savefig('SteadyState_dpDensity.png', dpi=1000)
plt.close()

for i in range(0, len(txtFiles)):
    # Monodisperse B, use activity
    if xA[i] == 0:
        plt.scatter(peB[i], ssGPD[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Monodisperse A, use activity
    elif xA[i] == 100:
        plt.scatter(peA[i], ssGPD[i], label=str(peA[i]), c='k')
        plt.xlabel(r'Activity')
    # Binary, varying xA
    elif xA[i] != xA[i-1]:
        plt.scatter(xA[i], ssGPD[i], c='k')
        plt.xlabel(r'Particle Fraction $x_{A}$')
        plt.xlim(0,1)
    # Binary, use activity ratio
    else:
        ratio = float(peA[i] / peB[i])
        plt.scatter(ratio, ssGPD[i], label=str(ratio), c='k')
        plt.xlabel(r'Activity Ratio')

plt.ylabel('Steady-State Gas Phase Density')
plt.savefig('SteadyState_gpDensity.png', dpi=1000)
plt.close()

max_stars_repo_name: kolbt/whingdingdilly

max_issues_repo_name: kolbt/whingdingdilly

max_forks_repo_name: kolbt/whingdingdilly