Remove unmapped reads from SAM output

Snakefile

@@ -7,21 +7,22 @@ from snakemake.utils import min_version
 ##################################
 
 
+extensions=["fa", "fasta", "fq", "fastq"]
+
 def multiglob(patterns):
     files = []
     for pattern in patterns:
         files.extend(glob.glob(pattern))
+    files = list(map(Path, files))
     return files
 
 
 def get_all_query_filepaths():
-    return multiglob(
-        ["queries/*.fa", "queries/*.fq", "queries/*.fasta", "queries/*.fastq"]
-    )
+    return multiglob(expand("queries/*.{ext}", ext=extensions))
 
 
 def get_all_query_filenames():
-    return (Path(file).with_suffix("").name for file in get_all_query_filepaths())
+    return (file.with_suffix("").name for file in get_all_query_filepaths())
 
 
 def get_filename_for_all_queries():
@@ -43,11 +44,11 @@ batches = [x.strip() for x in open(config["batches"])]
 print(f"Batches: {batches}")
 
 qfiles = get_all_query_filepaths()
-print(f"Query files: {qfiles}")
+print(f"Query files: {list(map(str, qfiles))}")
 
-assemblies_dir = f"{config['download_dir']}/asms"
-cobs_dir = f"{config['download_dir']}/cobs"
-decompression_dir = config.get("decompression_dir", "intermediate/00_cobs")
+assemblies_dir = Path(f"{config['download_dir']}/asms")
+cobs_dir = Path(f"{config['download_dir']}/cobs")
+decompression_dir = Path(config.get("decompression_dir", "intermediate/00_cobs"))
 
 
 wildcard_constraints:
@@ -155,13 +156,18 @@ rule download_cobs_batch:
 ##################################
 ## Processing rules
 ##################################
+def get_query_file(wildcards):
+    query_file = multiglob(expand(f"queries/{wildcards.qfile}.{{ext}}", ext=extensions))
+    assert len(query_file) == 1
+    return query_file[0]
+
 rule fix_query:
     """Fix query to expected COBS format: single line fastas composed of ACGT bases only
     """
     output:
         fixed_query="intermediate/fixed_queries/{qfile}.fa",
     input:
-        original_query="queries/{qfile}.fa",
+        original_query=get_query_file,
     threads: 1
     conda:
         "envs/seqtk.yaml"
@@ -179,7 +185,7 @@ rule concatenate_queries:
     """Concatenate all queries into a single file, so we just need to run COBS/minimap2 just once per batch
     """
     output:
-        concatenated_query=f"intermediate/fixed_queries/{get_filename_for_all_queries()}.fa",
+        concatenated_query=f"intermediate/concatenated_query/{get_filename_for_all_queries()}.fa",
     input:
         all_queries=expand(
             "intermediate/fixed_queries/{qfile}.fa", qfile=get_all_query_filenames()
@@ -216,7 +222,7 @@ rule run_cobs:
         match=protected("intermediate/01_match/{batch}____{qfile}.xz"),
     input:
         cobs_index=f"{decompression_dir}/{{batch}}.cobs_classic",
-        fa="intermediate/fixed_queries/{qfile}.fa",
+        fa="intermediate/concatenated_query/{qfile}.fa",
     threads: config["cobs_thr"]  # Small number in order to guarantee Snakemake parallelization
     params:
         kmer_thres=config["cobs_kmer_thres"],
@@ -244,7 +250,7 @@ rule decompress_and_run_cobs:
         match=protected("intermediate/01_match/{batch}____{qfile}.xz"),
     input:
         compressed_cobs_index=f"{cobs_dir}/{{batch}}.cobs_classic.xz",
-        fa="intermediate/fixed_queries/{qfile}.fa",
+        fa="intermediate/concatenated_query/{qfile}.fa",
     threads: config["cobs_thr"]  # Small number in order to guarantee Snakemake parallelization
     params:
         kmer_thres=config["cobs_kmer_thres"],
@@ -278,7 +284,7 @@ rule translate_matches:
     output:
         fa="intermediate/02_filter/{qfile}.fa",
     input:
-        fa="intermediate/fixed_queries/{qfile}.fa",
+        fa="intermediate/concatenated_query/{qfile}.fa",
         all_matches=[
             f"intermediate/01_match/{batch}____{{qfile}}.xz" for batch in batches
         ],