PR Part 1/3 of Cookiecutter pre-commit workflow

.flake8

@@ -0,0 +1,11 @@
+[flake8]
+exclude =
+    .git,
+    __pycache__,
+    build,
+    dist,
+    doc/source/conf.py
+max-line-length = 115
+# Ignore some style 'errors' produced while formatting by 'black'
+# https://black.readthedocs.io/en/stable/guides/using_black_with_other_tools.html#labels-why-pycodestyle-warnings
+extend-ignore = E203

diffpy/snmf/__init__.py

@@ -17,9 +17,9 @@
 """
 
 # obtain version information
-__version__ = '0.0.1'
+__version__ = "0.0.1"
 
 # top-level import
-#from diffpy.snmf.pdfmorph_api import pdfmorph, morph_default_config, plot_morph
+# from diffpy.snmf.pdfmorph_api import pdfmorph, morph_default_config, plot_morph
 
 # End of file

diffpy/snmf/containers.py

@@ -40,17 +40,21 @@ def apply_stretch(self, m):
           stretching operation, and one vector is the second derivative of the stretching operation.
         """
         normalized_grid = np.arange(len(self.grid))
-        func = lambda stretching_factor: np.interp(normalized_grid / stretching_factor, normalized_grid, self.iq,
-                                                   left=0, right=0)
+        func = lambda stretching_factor: np.interp(
+            normalized_grid / stretching_factor, normalized_grid, self.iq, left=0, right=0
+        )
         derivative_func = numdifftools.Derivative(func)
         second_derivative_func = numdifftools.Derivative(derivative_func)
 
         stretched_component = func(self.stretching_factors[m])
         stretched_component_gra = derivative_func(self.stretching_factors[m])
         stretched_component_hess = second_derivative_func(self.stretching_factors[m])
 
-        return np.asarray(stretched_component), np.asarray(stretched_component_gra), np.asarray(
-            stretched_component_hess)
+        return (
+            np.asarray(stretched_component),
+            np.asarray(stretched_component_gra),
+            np.asarray(stretched_component_hess),
+        )
 
     def apply_weight(self, m, stretched_component=None):
         """Applies as weight factor to a component signal.

diffpy/snmf/io.py

@@ -42,11 +42,21 @@ def initialize_variables(data_input, number_of_components, data_type, sparsity=1
     signal_length = data_input.shape[0]
     number_of_signals = data_input.shape[1]
 
-    diagonals = [np.ones(number_of_signals - 2), -2 * np.ones(number_of_signals - 2), np.ones(number_of_signals - 2)]
-    smoothness_term = .25 * scipy.sparse.diags(diagonals, [0, 1, 2], shape=(number_of_signals - 2, number_of_signals))
+    diagonals = [
+        np.ones(number_of_signals - 2),
+        -2 * np.ones(number_of_signals - 2),
+        np.ones(number_of_signals - 2),
+    ]
+    smoothness_term = 0.25 * scipy.sparse.diags(
+        diagonals, [0, 1, 2], shape=(number_of_signals - 2, number_of_signals)
+    )
 
     hessian_helper_matrix = scipy.sparse.block_diag([smoothness_term.T @ smoothness_term] * number_of_components)
-    sequence = np.arange(number_of_signals * number_of_components).reshape(number_of_components, number_of_signals).T.flatten()
+    sequence = (
+        np.arange(number_of_signals * number_of_components)
+        .reshape(number_of_components, number_of_signals)
+        .T.flatten()
+    )
     hessian_helper_matrix = hessian_helper_matrix[sequence, :][:, sequence]
 
     return {
@@ -57,7 +67,7 @@ def initialize_variables(data_input, number_of_components, data_type, sparsity=1
         "smoothness": smoothness,
         "sparsity": sparsity,
         "smoothness_term": smoothness_term,
-        "hessian_helper_matrix": hessian_helper_matrix
+        "hessian_helper_matrix": hessian_helper_matrix,
     }
 
 

diffpy/snmf/optimizers.py

@@ -20,7 +20,7 @@ def get_weights(stretched_component_gram_matrix, linear_coefficient, lower_bound
       Has length C.
 
     lower_bound: 1d array like
-      The lower bound on the values of the output weights. Has the same dimensions of the function output. Each 
+      The lower bound on the values of the output weights. Has the same dimensions of the function output. Each
       element in 'lower_bound' determines the minimum value the corresponding element in the function output may take.
 
     upper_bound: 1d array like
@@ -43,8 +43,9 @@ def get_weights(stretched_component_gram_matrix, linear_coefficient, lower_bound
     solution_variable = cvxpy.Variable(problem_size)
 
     objective = cvxpy.Minimize(
-        linear_coefficient.T @ solution_variable + 0.5 * cvxpy.quad_form(solution_variable,
-                                                                         stretched_component_gram_matrix))
+        linear_coefficient.T @ solution_variable
+        + 0.5 * cvxpy.quad_form(solution_variable, stretched_component_gram_matrix)
+    )
     constraints = [lower_bound <= solution_variable, solution_variable <= upper_bound]
 
     cvxpy.Problem(objective, constraints).solve()

diffpy/snmf/stretchednmfapp.py

@@ -5,25 +5,47 @@
 from diffpy.snmf.containers import ComponentSignal
 from diffpy.snmf.io import load_input_signals, initialize_variables
 
-ALLOWED_DATA_TYPES = ['powder_diffraction', 'pd', 'pair_distribution_function', 'pdf']
+ALLOWED_DATA_TYPES = ["powder_diffraction", "pd", "pair_distribution_function", "pdf"]
 
 
 def create_parser():
     parser = argparse.ArgumentParser(
-        prog="stretched_nmf",
-        description="Stretched Nonnegative Matrix Factorization"
+        prog="stretched_nmf", description="Stretched Nonnegative Matrix Factorization"
     )
-    parser.add_argument('-i', '--input-directory', type=str, default=None,
-                        help="Directory containing experimental data. Defaults to current working directory.")
-    parser.add_argument('-o', '--output-directory', type=str,
-                        help="The directory where the results will be written. Defaults to '<input_directory>/snmf_results'.")
-    parser.add_argument('t', '--data-type', type=str, default=None, choices=ALLOWED_DATA_TYPES,
-                        help="The type of the experimental data.")
-    parser.add_argument('-l', '--lift-factor', type=float, default=1,
-                        help="The lifting factor. Data will be lifted by lifted_data = data + abs(min(data) * lift). Default is 1.")
-    parser.add_argument('number-of-components', type=int,
-                        help="The number of component signals for the NMF decomposition. Must be an integer greater than 0")
-    parser.add_argument('-v', '--version', action='version', help='Print the software version number')
+    parser.add_argument(
+        "-i",
+        "--input-directory",
+        type=str,
+        default=None,
+        help="Directory containing experimental data. Defaults to current working directory.",
+    )
+    parser.add_argument(
+        "-o",
+        "--output-directory",
+        type=str,
+        help="The directory where the results will be written. Defaults to '<input_directory>/snmf_results'.",
+    )
+    parser.add_argument(
+        "t",
+        "--data-type",
+        type=str,
+        default=None,
+        choices=ALLOWED_DATA_TYPES,
+        help="The type of the experimental data.",
+    )
+    parser.add_argument(
+        "-l",
+        "--lift-factor",
+        type=float,
+        default=1,
+        help="The lifting factor. Data will be lifted by lifted_data = data + abs(min(data) * lift). Default is 1.",
+    )
+    parser.add_argument(
+        "number-of-components",
+        type=int,
+        help="The number of component signals for the NMF decomposition. Must be an integer greater than 0",
+    )
+    parser.add_argument("-v", "--version", action="version", help="Print the software version number")
     args = parser.parse_args()
     return args
 
@@ -35,5 +57,5 @@ def main():
     grid, input_data = load_input_signals(args.input_directory)
     lifted_input_data = lift_data(input_data, args.lift_factor)
     variables = initialize_variables(lifted_input_data, args.number_of_components, args.data_type)
-    components = initialize_components(variables['number_of_components'], variables['number_of_signals'], grid)
+    components = initialize_components(variables["number_of_components"], variables["number_of_signals"], grid)
     return components

diffpy/snmf/subroutines.py

@@ -170,11 +170,15 @@ def update_weights(components, data_input, method=None):
         stretched_components = np.zeros((signal_length, number_of_components))
         for i, component in enumerate(components):
             stretched_components[:, i] = component.apply_stretch(signal)[0]
-        if method == 'align':
+        if method == "align":
             weights = lsqnonneg(stretched_components, data_input[:, signal])
         else:
-            weights = get_weights(stretched_components.T @ stretched_components,
-                                  -stretched_components.T @ data_input[:, signal], 0, 1)
+            weights = get_weights(
+                stretched_components.T @ stretched_components,
+                -stretched_components.T @ data_input[:, signal],
+                0,
+                1,
+            )
             weight_matrix[:, signal] = weights
     return weight_matrix
 
@@ -236,13 +240,16 @@ def initialize_arrays(number_of_components, number_of_moments, signal_length):
     """
     component_matrix_guess = np.random.rand(signal_length, number_of_components)
     weight_matrix_guess = np.random.rand(number_of_components, number_of_moments)
-    stretching_matrix_guess = np.ones(number_of_components, number_of_moments) + np.random.randn(number_of_components,
-                                                                                                 number_of_moments) * 1e-3
+    stretching_matrix_guess = (
+        np.ones(number_of_components, number_of_moments)
+        + np.random.randn(number_of_components, number_of_moments) * 1e-3
+    )
     return component_matrix_guess, weight_matrix_guess, stretching_matrix_guess
 
 
-def objective_function(residual_matrix, stretching_factor_matrix, smoothness, smoothness_term, component_matrix,
-                       sparsity):
+def objective_function(
+    residual_matrix, stretching_factor_matrix, smoothness, smoothness_term, component_matrix, sparsity
+):
     """Defines the objective function of the algorithm and returns its value.
 
     Calculates the value of '(||residual_matrix||_F) ** 2 + smoothness * (||smoothness_term *
@@ -284,8 +291,11 @@ def objective_function(residual_matrix, stretching_factor_matrix, smoothness, sm
     residual_matrix = np.asarray(residual_matrix)
     stretching_factor_matrix = np.asarray(stretching_factor_matrix)
     component_matrix = np.asarray(component_matrix)
-    return .5 * np.linalg.norm(residual_matrix, 'fro') ** 2 + .5 * smoothness * np.linalg.norm(
-        smoothness_term @ stretching_factor_matrix.T, 'fro') ** 2 + sparsity * np.sum(np.sqrt(component_matrix))
+    return (
+        0.5 * np.linalg.norm(residual_matrix, "fro") ** 2
+        + 0.5 * smoothness * np.linalg.norm(smoothness_term @ stretching_factor_matrix.T, "fro") ** 2
+        + sparsity * np.sum(np.sqrt(component_matrix))
+    )
 
 
 def get_stretched_component(stretching_factor, component, signal_length):
@@ -325,11 +335,23 @@ def stretched_component_func(stretching_factor):
     stretched_component_gra = derivative_func(stretching_factor)
     stretched_component_hess = second_derivative_func(stretching_factor)
 
-    return np.asarray(stretched_component), np.asarray(stretched_component_gra), np.asarray(stretched_component_hess)
-
-
-def update_weights_matrix(component_amount, signal_length, stretching_factor_matrix, component_matrix, data_input,
-                          moment_amount, weights_matrix, method):
+    return (
+        np.asarray(stretched_component),
+        np.asarray(stretched_component_gra),
+        np.asarray(stretched_component_hess),
+    )
+
+
+def update_weights_matrix(
+    component_amount,
+    signal_length,
+    stretching_factor_matrix,
+    component_matrix,
+    data_input,
+    moment_amount,
+    weights_matrix,
+    method,
+):
     """Update the weight factors matrix.
 
     Parameters
@@ -376,21 +398,25 @@ def update_weights_matrix(component_amount, signal_length, stretching_factor_mat
     for i in range(moment_amount):
         stretched_components = np.zeros((signal_length, component_amount))
         for n in range(component_amount):
-            stretched_components[:, n] = get_stretched_component(stretching_factor_matrix[n, i], component_matrix[:, n],
-                                                                 signal_length)[0]
-        if method == 'align':
+            stretched_components[:, n] = get_stretched_component(
+                stretching_factor_matrix[n, i], component_matrix[:, n], signal_length
+            )[0]
+        if method == "align":
             weight = lsqnonneg(stretched_components[0:signal_length, :], data_input[0:signal_length, i])
         else:
             weight = get_weights(
                 stretched_components[0:signal_length, :].T @ stretched_components[0:signal_length, :],
                 -1 * stretched_components[0:signal_length, :].T @ data_input[0:signal_length, i],
-                0, 1)
+                0,
+                1,
+            )
         weights_matrix[:, i] = weight
     return weights_matrix
 
 
-def get_residual_matrix(component_matrix, weights_matrix, stretching_matrix, data_input, moment_amount,
-                        component_amount, signal_length):
+def get_residual_matrix(
+    component_matrix, weights_matrix, stretching_matrix, data_input, moment_amount, component_amount, signal_length
+):
     """Obtains the residual matrix between the experimental data and calculated data
 
     Calculates the difference between the experimental data and the reconstructed experimental data created from the
@@ -442,9 +468,11 @@ def get_residual_matrix(component_matrix, weights_matrix, stretching_matrix, dat
     for m in range(moment_amount):
         residual = residual_matrx[:, m]
         for k in range(component_amount):
-            residual = residual + weights_matrix[k, m] * get_stretched_component(stretching_matrix[k, m],
-                                                                                 component_matrix[:, k], signal_length)[
-                0]
+            residual = (
+                residual
+                + weights_matrix[k, m]
+                * get_stretched_component(stretching_matrix[k, m], component_matrix[:, k], signal_length)[0]
+            )
         residual_matrx[:, m] = residual
     return residual_matrx