begin snippy

pipeline.nf

@@ -35,8 +35,10 @@ def helpMessage() {
 // Extra configuration variables
 // SQLite commands script
 params.sqlite_commands = "$baseDir/sqlite_import.sh"
-params.genbank_asm_suffix = "_genomic.fna.gz"
+params.genbank_asm_gz_suffix = "_genomic.fna.gz"
+params.genbank_asm_fna_suffix = "_genomic.fna"
 params.assembly_for_download_ftp_file = "assembly_for_download.txt"
+params.reference_genome_ftp = "ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/009/065/GCF_000009065.1_ASM906v1/GCF_000009065.1_ASM906v1_genomic.fna.gz"
 
 // Show help message
 params.help = false
@@ -80,11 +82,11 @@ if (params.sqlite){
     """
     ftp_list=""
     fna_gz_list=""
-    echo "ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/009/909/655/GCA_009909655.1_ASM990965v1 ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/009/909/645/GCA_009909645.1_ASM990964v1" | sed 's/ /\\n/g' | while read line;
+    sqlite3 ${sqlite} 'select AssemblyFTPGenbank from Assembly' | head -n 50 | sed 's/ /\\n/g' | while read line;
     do
       if [[ ! -z \$line ]]; then
         asm_url=\$line;
-        asm_fasta=`echo \$line | cut -d "/" -f 10 | awk -v suffix=${params.genbank_asm_suffix} '{print \$0 suffix}'`;
+        asm_fasta=`echo \$line | cut -d "/" -f 10 | awk -v suffix=${params.genbank_asm_gz_suffix} '{print \$0 suffix}'`;
         asm_ftp=\${asm_url}/\${asm_fasta};
         echo \$asm_ftp >> ${params.assembly_for_download_ftp_file}
       fi;
@@ -97,22 +99,52 @@ if (params.sqlite){
     publishDir "${params.outdir}/assembly_download", mode: 'copy'
 
     echo true
-    // Deal with new lines
+
+    // Deal with new lines, split up ftp links by url
+    // By loading with file(), stages as local file
     assembly_for_download_ftp_ch.splitText()
             .map { file(it.replaceFirst(/\n/,'')) }
             .set { ftp_url_list_ch }
 
     input:
     file asm_fna_gz from ftp_url_list_ch
 
+    output:
+    file "*${params.genbank_asm_fna_suffix}" into asm_fna_ch
 
     script:
     """
-    echo ${asm_fna_gz}
+    # Use -f otherwise error due to too many levels of symbolic links
+    gunzip -f ${asm_fna_gz}
     """
   }
 }
 
+process snippy_pairwise{
+  // Pairwise align contigs to reference genome with snippy
+  echo true
+
+  input:
+  file asm_fna from asm_fna_ch
+  file reference_genome_fna from file(params.reference_genome_ftp)
+
+  output:
+
+  script:
+  """
+  snippy \
+    --cpus 2 \
+    --reference ${reference_genome_fna} \
+    --outdir output10X/${asm_fna.baseName} / \
+    --ctgs ${asm_fna} \
+    --mapqual 30 \
+    --mincov 10 \
+    --minfrac 0.9 \
+    --basequal 20 \
+    --report
+  """
+}
+
 
 
 def pipelineHeader() {