#!/usr/local/bin/perl
###############################################################################
# Program : ShowSWATHProphetChromatogram
#
# Description : This CGI program displays a single chromatogram in PeptideAtlas
#
#
###############################################################################
###############################################################################
# Basic SBEAMS setup
###############################################################################
use strict;
use FindBin;
use File::Basename;
use Data::Dumper;
use Carp;
use JSON;
$SIG{__DIE__} = sub { &Carp::confess };
use CGI::Carp qw (fatalsToBrowser);
use lib "$FindBin::Bin/../../lib/perl";
use vars qw ($q $sbeams $sbeamsMOD $PROG_NAME
$current_username $massCalculator );
use SBEAMS::Connection qw($q);
use SBEAMS::Connection::Settings;
use SBEAMS::Connection::Tables;
use SBEAMS::Connection::DataTable;
use SBEAMS::Proteomics::ChromatogramViewer;
use SBEAMS::Proteomics::PeptideMassCalculator;
use SBEAMS::PeptideAtlas;
use SBEAMS::PeptideAtlas::Settings;
use SBEAMS::PeptideAtlas::Tables;
use SBEAMS::PeptideAtlas::Chromatogram;
#$q = new CGI;
$sbeams = new SBEAMS::Connection;
$sbeamsMOD = new SBEAMS::PeptideAtlas;
$sbeamsMOD->setSBEAMS($sbeams);
$PROG_NAME="ShowChromatogram";
# Allowable neutral losses
my @neutral_loss_masses = ();
###############################################################################
# Define global variables if any and execute main()
###############################################################################
main();
###############################################################################
# Main Program:
#
# If $sbeams->Authenticate() succeeds, print header, process the CGI request,
# print the footer, and end.
###############################################################################
sub main {
#### Do the SBEAMS authentication and exit if a username is not returned
exit unless ($current_username = $sbeams->Authenticate(
permitted_work_groups_ref=>['PeptideAtlas_user','PeptideAtlas_admin',
'PeptideAtlas_readonly', 'PeptideAtlas_exec'],
#connect_read_only=>1,
allow_anonymous_access=>1,
));
#### Read in the default input parameters
my %parameters;
my $n_params_found =
$sbeams->parse_input_parameters(q=>$q,parameters_ref=>\%parameters);
#### Process generic "state" parameters before we start
$sbeams->processStandardParameters(parameters_ref=>\%parameters);
my $apply_action = $q->param('apply_action');
#### Process certain actions, then print the header, figure and do what the user wants, and print footer
$sbeamsMOD->display_page_header( show_doctype => 1 );
#processRequest();
processRequest(ref_parameters=>\%parameters);
$sbeamsMOD->display_page_footer();
$sbeams->display_page_footer(close_tables=>'YES',
separator_bar=>'YES',display_footer=>'NO');
} # end main
###############################################################################
# processRequest
###############################################################################
sub processRequest {
#### Process the arguments list
my %args = @_;
my $ref_parameters = $args{'ref_parameters'}
|| die "ref_parameters not passed";
$ref_parameters->{show_sptxt_info} = 1 if !defined $ref_parameters->{show_sptxt_info};
my $linkage = qq~
lgs[+87]VQAPSYGAR2
le[+87]AAIADAEQR2
anlet[-128.1]AIADAEQR2
avl[-71]EYLTAEILELAGNAAR3
VTIAQGGVLPNIqa[-71]2
LVl[-113.1]EVDPNIQAVR2
VTIAQGGVLPNIqav[-99.1]2:
VTIAQGGVLPNIq[+71]2
VTIAQGGVLPNIqa[+99.1]2
VTIAQGGVLPNIq[+170.1]2
gill[-57]HLESELAQTR2
EGNASGVSLLEALDt[+61]ILPPTRPTDKPLR4
lv[+97.1]LEVDPNIQAVR2
gill[-14]HLESELAQTR2
gill[+114]HLESELAQTR2
~;
#### Define param usage and print if parameter "help=1" defined.
my $usage = qq~
file - MRML file to pull from
pep - peptide sequence with odd mod nomenclature, with charge and something else appended (lgs[+87]VQAPSYGAR2:10)
rank - peak group rank (optional) 1
times - range of times 1391.8-1470.5
frags - fragments expected, with charge y8_1,y6_1,y7_1,b6_1,b4_1,b5_1
fragtimes - fragment elution times, for caption 1397.5-1453.7,1397.5-1464.9,1397.5-1453.7,1391.8-1456.5,1397.5-1470.5,1394.6-1470.5
prob - peak group probability (optional) 1
~;
$sbeams->show_help_if_requested(
usage_string => $usage,
ref_parameters => $ref_parameters,
);
my @frags = split( /,/, $ref_parameters->{frags} );
my @fragtimes = split( /,/, $ref_parameters->{fragtimes} );
my %frag_annot;
if ( @frags && @fragtimes ) {
for ( my $i = 0; $i <= $#frags; $i++ ) {
$frag_annot{$frags[$i]} = $fragtimes[$i];
}
}
# Create a chromatogram object so we can use its methods
my $cgram = new SBEAMS::PeptideAtlas::Chromatogram;
$cgram->setSBEAMS($sbeams);
# Fixme
my $extract_exe = '/net/dblocal/wwwSSL/html/devDC/sbeams/lib/scripts/PeptideAtlas/extractFromMRML.pl';
my $spectrum_pathname = $ref_parameters->{'spectrum_pathname'} || $ref_parameters->{file} || die "missing param spectrum pathname";
my $parent = 'PARENT';
my $seqstr = $ref_parameters->{'pep'};
my ( $seqchg, $sbin ) = split( /:/, $seqstr );
$seqchg =~ /^(.*)(\d)$/;
my $pepseq = $1;
my $chg = $2;
die "missing or improper parameter pep" if ( !$pepseq || !$chg || !$sbin );
my $time_range = $ref_parameters->{'times'};
# my $tdelta = abs(eval( "$time_range" ));
$time_range =~ s/-/,/;
my $command = "$extract_exe $spectrum_pathname $pepseq $chg RANGE=$time_range $parent";
my @results = `$command`;
my $rank = ( $ref_parameters->{rank} ) ? "Rank: $ref_parameters->{rank}" : '';
my $prob = ( $ref_parameters->{prob} ) ? "[Probability: $ref_parameters->{prob}]" : '';
print $linkage;
if ( !scalar( @results ) ) {
die "No results returned from given parameters
" . Dumper( $ref_parameters );
}
my @headings;
my %results;
my $extra;
my $rtable = SBEAMS::Connection::DataTable->new( BORDER => 1 );
for my $line ( @results ) {
chomp $line;
$line =~ s/\r//g;
chomp $line;
my @line = split( /\t/, $line, -1 );
$rtable->addRow( \@line ) unless $line =~ /PRED_INTS/;
if ( !scalar( @headings ) ) {
@headings = @line;
$headings[0] =~ s/^#//g;
next;
}
if ( $line =~ /PRED_INTS/ ) {
$extra = $line;
next;
}
for ( my $idx = 0; $idx <= $#line; $idx++ ) {
$results{$headings[$idx]} ||= [];
push @{$results{$headings[$idx]}}, $line[$idx];
}
}
my %maxint;
if ( $extra ) {
if ( $extra =~ /PRED_INTS:\s*(.*)$/ ) {
my $pred_ints = $1;
my @maxima = split( /,/, $pred_ints );
for my $pair ( @maxima ) {
my ( $ion, $max ) = split( /:/, $pair );
$maxint{$ion} = $max;
}
} else {
die "Regex failed for $extra";
}
#$VAR1 = '# MAXINT: 3518.00,5193.00 PRED_INTS: b6_1:990.6,y4_1:900.4,y7_1:3773.2,b9_1:674.5,b7_1:706.2,y5_1:1216.4';
}
my @ion_order;
if ( scalar( keys( %maxint ) ) ) {
@ion_order = sort { $maxint{$b} <=> $maxint{$a} } ( keys( %maxint ) );
unshift @ion_order, 'parent';
} else {
@ion_order = sort ( @headings );
}
my $spacer = ' ' x 36;
my $to_json = { info => [ { top_html => "$spacer$pepseq:$chg $prob $rank SWATH bin:$sbin
" } ],
data_json => [] };
my ( $start, $stop ) = split( /,/, $time_range );
my $median_in_minutes = ($start+$stop)/120;
my @time = @{$results{'time'}};
my $json = JSON->new();
for my $label ( @ion_order ) {
next if $label eq 'time';
my @time = @{$results{'time'}};
if ( $results{$label} ) {
} else {
next;
}
my @inten = @{$results{$label}};
my @data;
if ( $ref_parameters->{expand_time_range} ) {
# push @data, { 'time' => $median_in_minutes - $ref_parameters->{expand_time_range}, intensity => 0 };
# push @data, { 'time' => $time[0] - 0.01, intensity => 0 };
}
# my $time_in_minutes = ( $time[$i] ) ? $time[$i]/60 : 0;
# push @data, { 'time' => $time_in_minutes, intensity => $inten[$i] + 0 };
for ( my $i = 0; $i <= $#time; $i++ ) {
my $time_in_minutes = ( $time[$i] ) ? $time[$i]/60 : 0;
push @data, { 'time' => $time_in_minutes, intensity => $inten[$i] + 0 };
}
# if ( $ref_parameters->{expand_time_range} ) {
# push @data, { 'time' => $time[$#time] + 0.01, intensity => 0 };
# push @data, { 'time' => $median_in_minutes + $ref_parameters->{expand_time_range}, intensity => 0 };
# }
my $full_label = $label;
if ( scalar( keys( %maxint ) ) ) {
if ( $maxint{$label} ) {
$full_label .= "($maxint{$label})"
}
if ( $frag_annot{$label} ) {
$full_label .= " ET:$frag_annot{$label})"
}
}
if ( $label eq 'parent' ) {
unshift @{$to_json->{data_json}}, { label => $full_label, data => \@data, full => 'COUNT: 01 Q1:0.000 Q3:0.000' }
} else {
push @{$to_json->{data_json}}, { label => $full_label, data => \@data, full => 'COUNT: 01 Q1:0.000 Q3:0.000' }
}
}
my $json_string = $json->encode( $to_json );
# "info" : [
# {
# "top_html" : " (1315.675 Daltons)
Peptide: ILDETLYENAK[136]
Instrument: ABI SCIEX QTRAP 5500 "
# }
# ],
# /proteomics/sbader/forAndy/input/URINE/20130522_neat_1-1/UR_LIB/FULLMODS/TEST/REDO/20130522_neat_1-1_SW16_mrml.xml 'ASDC[CAM]GAGPIGFAGTVR[HvR]' 2 RANGE=4177.1,4257.2 PARENT | cut -f1
# https://db.systemsbiology.net/devDC/sbeams/cgi/PeptideAtlas/ShowSWATHProphetChromatogram?pepseq=ILDETLYENAK[136];use_pepname=1;limit_smoothing_options=1;default_smoothing_factor=5;spectrum_pathname=/proteomics/dcampbel/OpenSwath/test/OpenSWATH_Tutorial/split_napedro_L120420_010_SW-400AQUA_combined.chrom.idx.mzML;peptide=ILDETLYENAK;no_specfile=1;no_mquest=1;instrument_name=ABI%20SCIEX%20QTRAP%205500
my $sptxt_fulltext = '';
# Fetch some of the parameters into scalar variables and do some processing.
my $precursor_charge = $ref_parameters->{'precursor_charge'};
my $pepseq = $ref_parameters->{'pepseq'};
my $modified_pepseq = $ref_parameters->{'modified_pepseq'} || $pepseq;
my $spectrum_pathname = $ref_parameters->{'file'} ||
die 'ShowChromatogram; Need parameter spectrum_pathname';
my $spectrum_basename = $ref_parameters->{'spectrum_basename'};
$spectrum_pathname =~ /.*(\..*)/; my $filename_extension = $1;
if (! $spectrum_basename ) {
$spectrum_pathname =~ /.*\/(\S+?)\.${filename_extension}/;
$spectrum_basename = $1;
}
my $precursor_neutral_mass = $ref_parameters->{'precursor_neutral_mass'};
if ($ref_parameters->{'isotype'} =~ /heavy/i) {
$precursor_neutral_mass += $ref_parameters->{'isotype_delta_mass'};
}
my $machine = $ref_parameters->{'machine'}; #unused 12/06/11
#### Display the Chromatogram Viewer page.
# Print the HTML for the top of the chromatogram viewer page
#--------------------------------------------------
# my $top_html = $cgram->getTopHTMLforChromatogramViewer (
# param_href => $ref_parameters,
# seq => $modified_pepseq,
# precursor_neutral_mass => $precursor_neutral_mass,
# precursor_charge => $precursor_charge,
# spectrum_pathname => $spectrum_pathname,
# );
#--------------------------------------------------
my $json_string_no_newlines = $json_string;
$json_string_no_newlines =~ s/\n/ /g;
$json_string_no_newlines =~ s/'/"/g;
my $top_html = $cgram->getTopHTMLfromJson (
json_string => $json_string_no_newlines,
);
print $top_html;
my %smoothing;
if ( $ref_parameters->{default_smoothing_factor} ) {
$smoothing{default_smoothing_factor} = $ref_parameters->{default_smoothing_factor};
}
if ( $ref_parameters->{limit_smoothing_options} ) {
$smoothing{limit_smoothing_options} = $ref_parameters->{limit_smoothing_options};
}
if ( $ref_parameters->{expand_timeframe} ) {
$smoothing{expand_timeframe} = $ref_parameters->{expand_timeframe};
}
# Generate the Chromavis javascript chromatogram viewer
my $chromavis = new SBEAMS::Proteomics::ChromatogramViewer;
print $chromavis->generateChromatogram(
chromatogram_id => $ref_parameters->{'SEL_chromatogram_id'},
spectrum_pathname => $spectrum_pathname,
precursor_neutral_mass=> $precursor_neutral_mass,
precursor_charge => $precursor_charge,
seq => $modified_pepseq,
precursor_rt => $ref_parameters->{rt},
best_peak_group_rt => $ref_parameters->{Tr},
m_score => $ref_parameters->{m_score},
json_string => $json_string,
%smoothing
);
# Read the data from the JSON object into a hash of the
# format required by writeResultSet
my $cgram_href = $cgram->readJsonChromatogramIntoResultsetHash (
param_href => $ref_parameters,
json_string => $json_string,
#json_physical_pathname => $json_physical_pathname,
);
# Store chromatogram data as a recallable resultset
my $rs_set_name = "SETME";
$sbeams->writeResultSet(
resultset_file_ref=>\$rs_set_name,
resultset_ref=>$cgram_href,
file_prefix=>'chromgram_',
query_parameters_ref=>$ref_parameters,
);
# Temp, print ref_params;
my $ptable;
for my $p ( keys( %{$ref_parameters} ) ){
$ptable .= "$p => $ref_parameters->{$p}
\n";
}
print "$ptable
";
print "Data from MRML extraction
\n";
print "$extra
";
print "$rtable
";
# Print sptxt info, if provided.
my ($SRMProb) = ($sptxt_fulltext =~ /SRMProb=(\S*)\s*/);
if ( $ref_parameters->{show_sptxt_info} ) {
print "
SRM Probability: $SRMProb
\n" if $SRMProb;
print "$sptxt_fulltext
\n";
}
# Print the HTML that allows user to recall and download resultset
print $cgram->getBottomHTMLforChromatogramViewer(
param_href => $ref_parameters,
rs_set_name => $rs_set_name,
);
} # end processRequest