From 513db8b99b152068a57db5987293b4d741576c04 Mon Sep 17 00:00:00 2001
From: Matt Mallory <matt.mallory@alleninstitute.org>
Date: Wed, 19 Mar 2025 15:28:03 -0700
Subject: [PATCH 1/3] proj matrix edits

---
 morph_utils/ccf.py                            | 33 ++++---
 .../aggregate_single_cell_projection_csvs.py  | 92 +++++++++++++++++++
 .../projection_matrix_for_single_cell.py      | 11 ++-
 .../projection_matrix_from_swc_directory.py   | 86 +++++++++++++++--
 setup.cfg                                     |  1 +
 5 files changed, 195 insertions(+), 28 deletions(-)
 create mode 100644 morph_utils/executable_scripts/aggregate_single_cell_projection_csvs.py

diff --git a/morph_utils/ccf.py b/morph_utils/ccf.py
index 46c7523..702841f 100644
--- a/morph_utils/ccf.py
+++ b/morph_utils/ccf.py
@@ -289,7 +289,7 @@ def get_ccf_structure(voxel, name_map=None, annotation=None, coordinate_to_voxel
     return name_map[structure_id]
 
 def projection_matrix_for_swc(input_swc_file, mask_method = "tip_and_branch", 
-                              tip_count = False, annotation=None, 
+                              count_method = "node", annotation=None, 
                               annotation_path = None, volume_shape=(1320, 800, 1140),
                               resolution=10, node_type_list=[2],
                               resample_spacing=None):
@@ -304,7 +304,8 @@ def projection_matrix_for_swc(input_swc_file, mask_method = "tip_and_branch",
         'tip_and_branch' will return a projection matrix masking only structures with tip and branch nodes. If 'tip'
         will only look at structures with tip nodes. And last, if 'branch' will only look at structures with 
         branch nodes.
-        tip_count (bool): if True, will count number of tips instead of the length of axon
+        count_method (str): ['node','tip','branch']. When 'node', will measure axon length by multiplying by internode spacing.
+        Otherwise will return the count of tip or branch nodes in each structure
         annotation (array, optional): 3 dimensional ccf annotation array. Defaults to None.
         annotation_path (str, optional): path to nrrd file to use (optional). Defaults to None.
         volume_shape (tuple, optional): the size in voxels of the ccf atlas (annotation volume). Defaults to (1320, 800, 1140).
@@ -317,7 +318,7 @@ def projection_matrix_for_swc(input_swc_file, mask_method = "tip_and_branch",
         filename (str)
         
         specimen_projection_summary (dict): keys are strings of structures and values are the quantitiave projection
-        values. Either axon length, or number numbe of nodes depending on tip_count argument.
+        values. Either axon length, or number numbe of nodes depending on count_method argument.
 
     """
     
@@ -329,8 +330,11 @@ def projection_matrix_for_swc(input_swc_file, mask_method = "tip_and_branch",
                 print(f"WARNING: Annotation path provided does not exist, defaulting to 10um resolution, (1320,800, 1140) ccf.\n{annotation_path}")
                 annotation_path = None
         annotation = open_ccf_annotation(with_nrrd=True, annotation_path=annotation_path)
-        
-
+    
+    if count_method not in ['node','tip','branch']:
+        msg = f"count_method must be  'node','tip', or 'branch'. You passed in: {count_method}"        
+        raise ValueError(msg)
+    
     sg_df = load_structure_graph()
     name_map = NAME_MAP
     full_name_to_abbrev_dict = dict(zip(sg_df.name, sg_df.index))
@@ -379,15 +383,13 @@ def node_ider(morph,i):
 
     # determine ipsi/contra projections
     morph_df["ccf_structure_sided"] = morph_df.apply(lambda row: "ipsi_{}".format(row.ccf_structure) if row.z<z_midline else "contra_{}".format(row.ccf_structure), axis=1)
-
-
+    print(morph_df.head())
     # mask the morphology dataframe accordinagly
     if mask_method!="None":
             
         keep_structs = []
         for struct, struct_df in morph_df.groupby("ccf_structure_sided"):
             node_types_in_struct = struct_df.node_type.unique().tolist()
-        
             if (mask_method == 'tip') & ("tip" in node_types_in_struct):
                 keep_structs.append(struct)
                 
@@ -404,22 +406,25 @@ def node_ider(morph,i):
         
     else:
         morph_df_masked = morph_df
-
+    print(keep_structs)
+    print(morph_df_masked.head())
     # remove ventral targets and out of brain 
     ventral_targs = ["ipsi_{}".format(v) for v in vs_acronyms] + ["contra_{}".format(v) for v in vs_acronyms]
     targets_to_remove = ["ipsi_Out Of Cortex", "ipsi_root","contra_Out Of Cortex", "contra_root"] + ventral_targs
     morph_df_masked = morph_df_masked[~morph_df_masked['ccf_structure_sided'].isin(targets_to_remove)]
-
+    
+    print(morph_df_masked.node_type.value_counts())
+    print("count method",count_method)
     # accomodate tip counting instead of axon length 
-    if tip_count:
-        morph_df_masked = morph_df_masked[morph_df_masked['node_type']=='tip']
+    if count_method != 'node':
+        morph_df_masked = morph_df_masked[morph_df_masked['node_type']==count_method]
         spacing = 1
-    
+    print(morph_df_masked.head())
     # qunatify projections per structure 
     n_nodes_per_structure = morph_df_masked.ccf_structure_sided.value_counts()
     axon_length_per_structure = n_nodes_per_structure*spacing
     specimen_projection_summary = axon_length_per_structure.to_dict()
-    
+    print(specimen_projection_summary)
     return input_swc_file, specimen_projection_summary   
 
  
diff --git a/morph_utils/executable_scripts/aggregate_single_cell_projection_csvs.py b/morph_utils/executable_scripts/aggregate_single_cell_projection_csvs.py
new file mode 100644
index 0000000..944670b
--- /dev/null
+++ b/morph_utils/executable_scripts/aggregate_single_cell_projection_csvs.py
@@ -0,0 +1,92 @@
+import os
+from tqdm import tqdm
+import pandas as pd
+import argschema as ags
+from morph_utils.ccf import projection_matrix_for_swc
+
+class IO_Schema(ags.ArgSchema):
+    output_directory = ags.fields.OutputDir(description="output directory")
+    output_projection_csv = ags.fields.OutputFile(description="output projection csv")
+    mask_method = ags.fields.Str(default="tip_and_branch",description = " 'tip_and_branch', 'branch', 'tip', or 'tip_or_branch' ")
+    projection_threshold = ags.fields.Int(default=0)
+    normalize_proj_mat = ags.fields.Boolean(default=True)
+
+
+def normalize_projection_columns_per_cell(input_df, projection_column_identifiers=['ipsi', 'contra']):
+    """
+    :param input_df:  input projection df
+    :param projection_column_identifiers: list of identifiers for projection columns. i.e. strings that identify projection columns from metadata columns
+    :return: normalized projection matrix
+    """
+    proj_cols = [c for c in input_df.columns if any([ider in c for ider in projection_column_identifiers])]
+    input_df[proj_cols] = input_df[proj_cols].fillna(0)
+
+    res = input_df[proj_cols].T / input_df[proj_cols].sum(axis=1)
+    input_df[proj_cols] = res.T
+
+    return input_df
+
+
+def main(output_directory,
+         output_projection_csv, 
+         projection_threshold,
+         mask_method,
+         normalize_proj_mat,
+         **kwargs):
+    
+    files_of_interest = [f for f in os.listdir(output_directory) if (f.endswith(".csv") and not f.endswith("_norm.csv")) ]
+    output_projection_csv = output_projection_csv.replace(".csv", f"_{mask_method}.csv")
+    
+    projection_records = {}
+    # branch_and_tip_projection_records = {}
+    for fn in files_of_interest:
+        df = pd.read_csv(os.path.join(output_directory, fn),index_col=0)
+        src_file = df.index[0]
+        fn = os.path.abspath(src_file)
+        
+        proj_records = df.loc[src_file].to_dict()
+        # brnch_tip_records = res[1]
+
+        projection_records[fn] = proj_records
+        # branch_and_tip_projection_records[fn] = brnch_tip_records
+
+    proj_df = pd.DataFrame(projection_records).T.fillna(0)
+    # proj_df_mask = pd.DataFrame(branch_and_tip_projection_records).T.fillna(0)
+
+    proj_df.to_csv(output_projection_csv)
+    # proj_df_mask.to_csv(output_projection_csv_tip_branch_mask)
+
+    if projection_threshold != 0:
+        output_projection_csv = output_projection_csv.replace(".csv",
+                                                              "{}thresh.csv".format(projection_threshold))
+        # output_projection_csv_tip_branch_mask = output_projection_csv_tip_branch_mask.replace(".csv",
+        #                                                                                       "{}thresh.csv".format(
+        #                                                                                           projection_threshold))
+
+        proj_df_arr = proj_df.values
+        proj_df_arr[proj_df_arr < projection_threshold] = 0
+        proj_df = pd.DataFrame(proj_df_arr, columns=proj_df.columns, index=proj_df.index)
+        proj_df.to_csv(output_projection_csv)
+
+        # proj_df_mask_arr = proj_df_mask.values
+        # proj_df_mask_arr[proj_df_mask_arr < projection_threshold] = 0
+        # proj_df_mask = pd.DataFrame(proj_df_mask_arr, columns=proj_df_mask.columns, index=proj_df_mask.index)
+        # proj_df_mask.to_csv(output_projection_csv_tip_branch_mask)
+
+    if normalize_proj_mat:
+        output_projection_csv = output_projection_csv.replace(".csv", "_norm.csv")
+        # output_projection_csv_tip_branch_mask = output_projection_csv_tip_branch_mask.replace(".csv", "_norm.csv")
+
+        proj_df = normalize_projection_columns_per_cell(proj_df)
+        proj_df.to_csv(output_projection_csv)
+
+        # proj_df_mask = normalize_projection_columns_per_cell(proj_df_mask)
+        # proj_df_mask.to_csv(output_projection_csv_tip_branch_mask)
+
+def console_script():
+    module = ags.ArgSchemaParser(schema_type=IO_Schema)
+    main(**module.args)
+
+if __name__ == "__main__":
+    module = ags.ArgSchemaParser(schema_type=IO_Schema)
+    main(**module.args)
diff --git a/morph_utils/executable_scripts/projection_matrix_for_single_cell.py b/morph_utils/executable_scripts/projection_matrix_for_single_cell.py
index ae1a844..057506d 100644
--- a/morph_utils/executable_scripts/projection_matrix_for_single_cell.py
+++ b/morph_utils/executable_scripts/projection_matrix_for_single_cell.py
@@ -10,7 +10,7 @@ class IO_Schema(ags.ArgSchema):
     projection_threshold = ags.fields.Int(default=0)
     normalize_proj_mat = ags.fields.Boolean(default=True)
     mask_method = ags.fields.Str(default="tip_and_branch",description = " 'tip_and_branch', 'branch', 'tip', or 'tip_or_branch' ")
-    tip_count = ags.fields.Boolean(default=False, description="when true, this will measure a matrix of number of tips instead of number of nodes")
+    count_method = ags.fields.String(default="node", description="should be a member of ['node','tip','branch']")
     annotation_path = ags.fields.Str(default="",description = "Optional. Path to annotation .nrrd file. Defaults to 10um ccf atlas")
     resolution = ags.fields.Int(default=10, description="Optional. ccf resolution (micron/pixel")
     volume_shape = ags.fields.List(ags.fields.Int, default=[1320, 800, 1140], description = "Optional. Size of input annotation")
@@ -38,7 +38,7 @@ def main(input_swc_file,
          projection_threshold, 
          normalize_proj_mat,
          mask_method,
-         tip_count,
+         count_method,
          annotation_path,
          volume_shape,
          resample_spacing,
@@ -46,10 +46,13 @@ def main(input_swc_file,
     
     if annotation_path == "":
         annotation_path = None
+        
+    if mask_method not in [None,'tip_and_branch', 'branch', 'tip', 'tip_or_branch']:
+        raise ValueError(f"Invalid mask_method provided {mask_method}")  
     
     results = []
     res = projection_matrix_for_swc(input_swc_file=input_swc_file, 
-                                    tip_count = tip_count, 
+                                    count_method = count_method, 
                                     mask_method = mask_method,
                                     annotation=None, 
                                     annotation_path = annotation_path, 
@@ -74,7 +77,7 @@ def main(input_swc_file,
 
     proj_df.to_csv(output_projection_csv)
     # proj_df_mask.to_csv(output_projection_csv_tip_branch_mask)
-
+    print(proj_df.head())
     if projection_threshold != 0:
         output_projection_csv = output_projection_csv.replace(".csv",
                                                               "{}thresh.csv".format(projection_threshold))
diff --git a/morph_utils/executable_scripts/projection_matrix_from_swc_directory.py b/morph_utils/executable_scripts/projection_matrix_from_swc_directory.py
index 8992f40..1c9e4b1 100644
--- a/morph_utils/executable_scripts/projection_matrix_from_swc_directory.py
+++ b/morph_utils/executable_scripts/projection_matrix_from_swc_directory.py
@@ -13,14 +13,14 @@ class IO_Schema(ags.ArgSchema):
     projection_threshold = ags.fields.Int(default=0)
     normalize_proj_mat = ags.fields.Boolean(default=True)
     mask_method = ags.fields.Str(default="tip_and_branch",description = " 'tip_and_branch', 'branch', 'tip', or 'tip_or_branch' ")
-    tip_count = ags.fields.Boolean(default=False, description="when true, this will measure a matrix of number of tips instead of number of nodes")
+    count_method = ags.fields.String(default="node", description="should be a member of ['node','tip','branch']")
     annotation_path = ags.fields.Str(default="",description = "Optional. Path to annotation .nrrd file. Defaults to 10um ccf atlas")
     resolution = ags.fields.Int(default=10, description="Optional. ccf resolution (micron/pixel")
     volume_shape = ags.fields.List(ags.fields.Int, default=[1320, 800, 1140], description = "Optional. Size of input annotation")
     resample_spacing = ags.fields.Float(allow_none=True, default=None, description = 'internode spacing to resample input morphology with')
     run_host = ags.fields.String(default='local',description='either ["local" or "hpc"]. Will run either locally or submit jobs to the HPC')
     virtual_env_name = ags.fields.Str(default='skeleton_keys_4',description='Name of virtual conda env to activate on hpc. not needed if running local')
-    output_projection_csv = ags.fields.OutputFile(allow_none=True,default=None, description="output projection csv, when running local only")
+    output_projection_csv = ags.fields.OutputFile(description="output projection csv, when running local only")
 
     
 def normalize_projection_columns_per_cell(input_df, projection_column_identifiers=['ipsi', 'contra']):
@@ -44,7 +44,7 @@ def main(ccf_swc_directory,
          projection_threshold, 
          normalize_proj_mat,
          mask_method,
-         tip_count,
+         count_method,
          annotation_path,
          volume_shape,
          resample_spacing,
@@ -55,6 +55,8 @@ def main(ccf_swc_directory,
     
     if run_host not in ['local','hpc']:
         raise ValueError(f"Invalid run_host parameter entered ({run_host})")
+    if mask_method not in [None,'tip_and_branch', 'branch', 'tip', 'tip_or_branch']:
+        raise ValueError(f"Invalid mask_method provided {mask_method}")
     
     if annotation_path == "":
         annotation_path = None
@@ -69,13 +71,14 @@ def main(ccf_swc_directory,
                 os.mkdir(dd)
                 
     results = []
+    single_cell_job_ids = []
     for swc_fn in tqdm([f for f in os.listdir(ccf_swc_directory) if ".swc" in f]):
 
         swc_pth = os.path.abspath(os.path.join(ccf_swc_directory, swc_fn))
         
         if run_host=='local':
             res = projection_matrix_for_swc(input_swc_file=swc_pth, 
-                                tip_count = tip_count,
+                                count_method= count_method,
                                 mask_method = mask_method,
                                 annotation=None, 
                                 annotation_path = annotation_path, 
@@ -86,24 +89,25 @@ def main(ccf_swc_directory,
 
         else:
                 
-            output_projection_csv = os.path.join(single_sp_proj_dir, swc_fn.replace(".swc",".csv"))
+            this_output_projection_csv = os.path.join(single_sp_proj_dir, swc_fn.replace(".swc",".csv"))
             
-            if not os.path.exists(output_projection_csv):
+            if not os.path.exists(this_output_projection_csv):
                     
                 job_file = os.path.join(job_dir,swc_fn.replace(".swc",".sh"))
                 log_file = os.path.join(job_dir,swc_fn.replace(".swc",".log"))
                 
                 command = "morph_utils_extract_projection_matrix_single_cell "
                 command = command+ f" --input_swc_file '{swc_pth}'"
-                command = command+ f" --output_projection_csv {output_projection_csv}"
+                command = command+ f" --output_projection_csv {this_output_projection_csv}"
                 command = command+ f" --projection_threshold {projection_threshold}"
                 command = command+ f" --normalize_proj_mat {normalize_proj_mat}"
                 command = command+ f" --mask_method {mask_method}"
-                command = command+ f" --tip_count {tip_count}"
+                command = command+ f" --count_method {count_method}"
                 command = command+ f" --annotation_path {annotation_path}"
                 command = command+ f" --resolution {resolution}"
                 # command = command+ f" --volume_shape {volume_shape}"
-                command = command+ f" --resample_spacing {resample_spacing}"
+                if resample_spacing is not None:
+                    command = command+ f" --resample_spacing {resample_spacing}"
 
                     
                 activate_command = f"source activate {virtual_env_name}"
@@ -150,8 +154,70 @@ def main(ccf_swc_directory,
                 std_out = result.stdout.decode('utf-8')
 
                 job_id = std_out.split("Submitted batch job ")[-1].replace("\n", "")
-                time.sleep(0.1)
+                single_cell_job_ids.append(job_id)
+                # time.sleep(0.1)
+                
+    if run_host!='local':
+        # aggregate single projection files into proj mat
+        job_file = os.path.join(job_dir,"Projection_Aggregation.sh")
+        log_file = os.path.join(job_dir,"Projection_Aggregation.log")
+        agg_outdir = os.path.join(output_directory,'SingleCellProjections')
+        command = "morph_utils_aggregate_single_cell_projs "
+        command = command+ f" --output_projection_csv {output_projection_csv}"
+        command = command+ f" --output_directory {agg_outdir}"
+        command = command+ f" --projection_threshold {projection_threshold}"
+        command = command+ f" --normalize_proj_mat {normalize_proj_mat}"
+        command = command+ f" --mask_method {mask_method}"
+            
+        activate_command = f"source activate {virtual_env_name}"
+        command_list = [activate_command, command]
+                
+        slurm_kwargs = {
+            "--job-name": "AggregateProjs",
+            "--mail-type": "NONE",
+            "--cpus-per-task": "1",
+            "--nodes": "1",
+            "--kill-on-invalid-dep": "yes",
+            "--mem": "4gb",
+            "--time": "1:00:00",
+            "--partition": "celltypes",
+            "--output": log_file
+        }
+
+        dag_node = {
+            "job_file":job_file,
+            "slurm_kwargs":slurm_kwargs,
+            "slurm_commands":command_list
+        }
+        
+        job_file = dag_node["job_file"]
+        slurm_kwargs = dag_node["slurm_kwargs"]
+        command_list = dag_node["slurm_commands"]
+
+        job_string_list = [f"#SBATCH {k}={v}" for k, v in slurm_kwargs.items()]
+        job_string_list = job_string_list + command_list
+        job_string_list = ["#!/bin/bash"] + job_string_list
+
+        if os.path.exists(job_file):
+            os.remove(job_file)
+
+        with open(job_file, 'w') as job_f:
+            for val in job_string_list:
+                job_f.write(val)
+                job_f.write('\n')
+                
+        command = "sbatch --dependency=afterany"
+        for p_jid in single_cell_job_ids:
+            command = command + f":{p_jid}"
+        command = command + " {}".format(job_file)
+        command_list = command.split(" ")
+        # print(command)
+        result = subprocess.run(command_list, stdout=subprocess.PIPE)
+        std_out = result.stdout.decode('utf-8')
+
+        job_id = std_out.split("Submitted batch job ")[-1].replace("\n", "")
         
+       
     if results != []:
         
         output_projection_csv = output_projection_csv.replace(".csv", f"_{mask_method}.csv")
diff --git a/setup.cfg b/setup.cfg
index 9b73ca8..c56fcad 100644
--- a/setup.cfg
+++ b/setup.cfg
@@ -26,5 +26,6 @@ console_scripts =
     morph_utils_full_morph_soma_correction = morph_utils.executable_scripts.full_morphology_soma_correction:console_script
     morph_utils_extract_projection_matrix = morph_utils.executable_scripts.projection_matrix_from_swc_directory:console_script
     morph_utils_extract_projection_matrix_single_cell = morph_utils.executable_scripts.projection_matrix_for_single_cell:console_script
+    morph_utils_aggregate_single_cell_projs = morph_utils.executable_scripts.aggregate_single_cell_projection_csvs:console_script
     morph_utils_move_somas_left_hemisphere = morph_utils.executable_scripts.move_somas_to_left_hemisphere_swc_directory:console_script
     morph_utils_local_crop_ccf_swcs = morph_utils.executable_scripts.local_crop_ccf_swc_directory:console_script

From 7028163b7d80b10328b235bf3e63286e00659cfc Mon Sep 17 00:00:00 2001
From: Matt Mallory <matt.mallory@alleninstitute.org>
Date: Mon, 19 May 2025 13:45:48 -0700
Subject: [PATCH 2/3] fix: small bug in resample

---
 morph_utils/modifications.py | 34 +++++++++++++++++++++++++---------
 1 file changed, 25 insertions(+), 9 deletions(-)

diff --git a/morph_utils/modifications.py b/morph_utils/modifications.py
index 01e7125..1773e32 100644
--- a/morph_utils/modifications.py
+++ b/morph_utils/modifications.py
@@ -129,9 +129,8 @@ def resample_morphology(morph, spacing_size):
         old_irr_id_to_new_irr_id_dict[root['id']]=new_root['id']
         this_roots_children = morph.get_children(root)
         node_ct+=1
-          
-        for child in this_roots_children:
-                
+        for child in this_roots_children:   
+            # in a sense children of the root nodes are also treated as irreducible nodes
             new_child = copy(child)
             new_child['parent'] = new_root['id']
             new_child['id'] = node_ct
@@ -151,18 +150,31 @@ def resample_morphology(morph, spacing_size):
                 
                 parent = morph.node_by_id(current_node['parent'])
                 siblings = morph.get_children(parent)
+                # so if the current_node has siblings (it's parent furcates),
+                # and we've already visited the parent, we should add that parent to this_list
+                # as it will be the first upstream irreducible node in this segment
                 if len(siblings)>1 and parent['id'] in seen_ids:
                     if parent['id']!=morph.get_soma()['id']:
+                        # this needs to happend before we add current_node to the list
                         this_list.append(parent)
-                    
+
+                # now add current node, and update that we've seen curent node
                 this_list.append(current_node)
                 seen_ids.update([current_node['id']])
                 children_list = morph.get_children(current_node)
                 if len(children_list)!=1:
-                    irreducible_segments.append(this_list)
+                    # if `current_node` is the immediate child of `child`, and is irreducible, 
+                    # then `this_list` will only contain `current_node`
+                    if this_list!= [current_node]:
+                        # otherwise, add this_list to our meta list
+                        irreducible_segments.append(this_list)
+                        
+                    # refresh the list
                     this_list = []
                 
                 for ch_no in children_list:
+                    # add the children, and the upstream irreducible parent node will be 
+                    # added using the logic from beore
                     queue.appendleft(ch_no)
 
 
@@ -227,8 +239,12 @@ def resample_morphology(morph, spacing_size):
                 node_ct+=1
                 new_node_2 = copy(irr_node_2)
                 if len(reducible_arr)==0:
+                    # if there are no reducible nodes, the down tree reducible node
+                    # should point to the upstream one
                     new_node_2['parent'] = red_1_parent_id
                 else:
+                    # otherwise, the downstream irreducible node should point to 
+                    # the last redubile node id, which will be node_ct-1
                     new_node_2['parent'] = node_ct-1
                 new_node_2['id']=node_ct
                 new_nodes.append(new_node_2)
@@ -239,10 +255,10 @@ def resample_morphology(morph, spacing_size):
 
             node_ct+=1    
 
-        resampled_morph = Morphology(new_nodes,
-                parent_id_cb=lambda x:x['parent'],
-                node_id_cb=lambda x:x['id'])
-        
+    resampled_morph = Morphology(new_nodes,
+            parent_id_cb=lambda x:x['parent'],
+            node_id_cb=lambda x:x['id'])
+
     return resampled_morph
 
 

From 6933b1227c367383f7d3e7de183b421e29f831eb Mon Sep 17 00:00:00 2001
From: Matt Mallory <matt.mallory@alleninstitute.org>
Date: Mon, 19 May 2025 13:53:53 -0700
Subject: [PATCH 3/3] remove prints

---
 morph_utils/ccf.py | 13 +++++--------
 1 file changed, 5 insertions(+), 8 deletions(-)

diff --git a/morph_utils/ccf.py b/morph_utils/ccf.py
index 702841f..79ea5dd 100644
--- a/morph_utils/ccf.py
+++ b/morph_utils/ccf.py
@@ -362,7 +362,9 @@ def projection_matrix_for_swc(input_swc_file, mask_method = "tip_and_branch",
     
     # annotate each node
     if morph_df.empty:
-        print("Its empty")
+        
+        msg = "morph_df is empty, possibly caused by `morph_df = morph_df[morph_df['type'].isin(node_type_list)]`"
+        warnings.warn(msg)
         return input_swc_file, {} 
     
     morph_df['ccf_structure'] = morph_df.apply(lambda rw: full_name_to_abbrev_dict[get_ccf_structure( np.array([rw.x, rw.y, rw.z]) , name_map, annotation, True)], axis=1)
@@ -383,7 +385,6 @@ def node_ider(morph,i):
 
     # determine ipsi/contra projections
     morph_df["ccf_structure_sided"] = morph_df.apply(lambda row: "ipsi_{}".format(row.ccf_structure) if row.z<z_midline else "contra_{}".format(row.ccf_structure), axis=1)
-    print(morph_df.head())
     # mask the morphology dataframe accordinagly
     if mask_method!="None":
             
@@ -405,26 +406,22 @@ def node_ider(morph,i):
         morph_df_masked = morph_df[morph_df['ccf_structure_sided'].isin(keep_structs)]
         
     else:
+        print("Not masking projection matrix...")
         morph_df_masked = morph_df
-    print(keep_structs)
-    print(morph_df_masked.head())
+        
     # remove ventral targets and out of brain 
     ventral_targs = ["ipsi_{}".format(v) for v in vs_acronyms] + ["contra_{}".format(v) for v in vs_acronyms]
     targets_to_remove = ["ipsi_Out Of Cortex", "ipsi_root","contra_Out Of Cortex", "contra_root"] + ventral_targs
     morph_df_masked = morph_df_masked[~morph_df_masked['ccf_structure_sided'].isin(targets_to_remove)]
     
-    print(morph_df_masked.node_type.value_counts())
-    print("count method",count_method)
     # accomodate tip counting instead of axon length 
     if count_method != 'node':
         morph_df_masked = morph_df_masked[morph_df_masked['node_type']==count_method]
         spacing = 1
-    print(morph_df_masked.head())
     # qunatify projections per structure 
     n_nodes_per_structure = morph_df_masked.ccf_structure_sided.value_counts()
     axon_length_per_structure = n_nodes_per_structure*spacing
     specimen_projection_summary = axon_length_per_structure.to_dict()
-    print(specimen_projection_summary)
     return input_swc_file, specimen_projection_summary