Skip to content

Commit

Permalink
analyze.php: Improved read data fallback to HG19 (now works if the ou…
Browse files Browse the repository at this point in the history
…tput folder exists as well)
  • Loading branch information
marc-sturm committed Dec 6, 2021
1 parent f026c65 commit 0a10a52
Show file tree
Hide file tree
Showing 2 changed files with 27 additions and 14 deletions.
37 changes: 25 additions & 12 deletions src/Pipelines/analyze.php
Original file line number Diff line number Diff line change
Expand Up @@ -31,10 +31,8 @@
extract($parser->parse($argv));

// create logfile in output folder if no filepath is provided:
$output_folder_exists = true;
if (!file_exists($folder))
{
$output_folder_exists = false;
exec2("mkdir -p $folder");
}
if ($parser->getLogFile() == "") $parser->setLogFile($folder."/analyze_".date("YmdHis").".log");
Expand Down Expand Up @@ -180,19 +178,34 @@
$bamfile_to_convert = $bamfile;

//fallback for remapping GRCh37 samples to GRCh38 - use BAM to create FASTQs
if (!$output_folder_exists && $sys['build']=="GRCh38")
if ($sys['build']=="GRCh38")
{
if (!db_is_enabled("NGSD")) trigger_error("NGSD access required to determine GRCh37 sample path!", E_USER_ERROR);
$db = DB::getInstance("NGSD", false);
$info = get_processed_sample_info($db, $name, false);
//determine if HG38 folder contains read data
$read_data_present = false;
list($files) = exec2("ls {$folder}");
foreach($files as $file)
{
if (ends_with($file, ".bam") || ends_with($file, ".fastq.gz"))
{
$read_data_present = true;
}
}

//no read data > try to generate it from HG19
if (!$read_data_present)
{
if (!db_is_enabled("NGSD")) trigger_error("NGSD access required to determine GRCh37 sample path!", E_USER_ERROR);
$db = DB::getInstance("NGSD", false);
$info = get_processed_sample_info($db, $name, false);

// determine project folder of GRCh37
$bamfile_to_convert = get_path("GRCh37_project_folder").$info['project_type']."/".$info['project_name']."/Sample_${name}/${name}.bam";
// determine project folder of GRCh37
$bamfile_to_convert = get_path("GRCh37_project_folder").$info['project_type']."/".$info['project_name']."/Sample_${name}/${name}.bam";

// check if bam file exists
if (!file_exists($bamfile_to_convert)) trigger_error("BAM file of GRCh37 sample is missing!", E_USER_ERROR);
// check if bam file exists
if (!file_exists($bamfile_to_convert)) trigger_error("BAM file of GRCh37 sample is missing!", E_USER_ERROR);

trigger_error("Output folder does not exists. Using BAM or Fastq files from GRCh37 as input!", E_USER_NOTICE);
trigger_error("No read data found in output folder. Using BAM/FASTQ files from GRCh37 as input!", E_USER_NOTICE);
}
}

//determine input FASTQ files
Expand Down Expand Up @@ -234,7 +247,7 @@
{
if(file_exists($bamfile_to_convert))
{
trigger_error("No FASTQ files found in folder. Using BAM file to generate FASTQ files.", E_USER_NOTICE);
trigger_error("No FASTQ files found in folder. Using BAM file to generate FASTQ files: $bamfile_to_convert", E_USER_NOTICE);

// extract reads from BAM file
$in_fq_for = $folder."/{$name}_BamToFastq_R1_001.fastq.gz";
Expand Down
4 changes: 2 additions & 2 deletions src/Tools/liftover_bed.php
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@
require_once(dirname($_SERVER['SCRIPT_FILENAME'])."/../Common/all.php");

// parse command line arguments
$parser = new ToolBase("liftover_bed", "Lift-over from.");
$parser = new ToolBase("liftover_bed", "Lift-over ob BED file from GRCh37 to GRCh38.");
$parser->addInfile("in", "Input BED file.", false);
$parser->addOutfile("out", "Input BED file.", false);
$parser->addInt("max_size_inc", "Maximum size increase.", true, 10);
Expand Down Expand Up @@ -54,7 +54,7 @@
$size_diff_perc = 100.0 * $size_diff / $size_before;
if ($size_after>$size_before && $size_diff>$max_size_inc && $size_diff_perc>$max_size_inc_perc)
{
print "Error - different size after mapping (from ".($end-$start)." to ".($end2-$start2)."): $line\n";
print "Error - bigger size after mapping (from ".($end-$start)." to ".($end2-$start2)."): $line\n";
continue;
}

Expand Down

0 comments on commit 0a10a52

Please sign in to comment.