package LCMerge;
use strict;
use English;
use Carp;
use List::Util qw(reduce);
use List::MoreUtils qw(uniq);
use vars
qw(@ISA $VERSION $name $author $date $version @instList $numberOfStreams);
@ISA = qw(Module);
use ModuleResources;
$name = __PACKAGE__;
$VERSION = '0.08';
$version = $VERSION;
$author = 'Jose Vicente Perea';
$date = '2017-06-06';
@instList = qw(rgs);
sub numberOfStreams
{
return 1;
}
sub evaluateRules
{
return ignore()
if allIgnored( module => 'EPICSourceProducts' ) and allIgnored( module => 'RGSProducts' );
start()
if (allComplete( module => 'EPICSourceProducts' ) and allComplete( module => 'RGSProducts' ));
}
sub performAction
{
info("Module version number: $version");
my $order = 1;
my %collect = ();
my (@rgsTSList, @rgs1TSList, @rgs2TSList, @epnTSList, @emos1TSList, @emos2TSList);
my (%rgs1, %rgs2);
my %src_files;
info("******************************************************************************************");
info("********** DEBUG - LCMerge - EPIC SRC_NUM identification in RGS Source Lists *************");
info("******************************************************************************************");
foreach my $inst ( "rgs2", "rgs1" )
{
foreach my $exp_id ( listExposures( instrument => $inst ) )
{
my $srcList = findFile( class => 'product'
, instrument => $inst
, exp_id => $exp_id
, content => 'RGS source list'
);
if (!defined($srcList) || $srcList eq '' || !fileExists(file => $srcList))
{
info("Can't find $inst source list for expid $exp_id");
next;
}
my %src_param = &ModuleUtil::get_rgssources_param( $srcList , 'SRCLIST' , $inst );
info("DEBUG - LCMerge - Inst: $inst , srcList: $srcList");
my $epic_src_num;
my @rgs_processed_sources = (keys %src_param);
info("DEBUG - LCMerge - Sources with RGS products, src_num_INDEX's : @rgs_processed_sources");
foreach my $src ( sort ( @rgs_processed_sources ) )
{
my $src_num = $src_param{$src}{INDEX};
my $src_label = $src_param{$src}{LABEL};
my $src_process = $src_param{$src}{PROCESS};
info("DEBUG - LCMerge -......");
info("DEBUG - LCMerge -......expid: $exp_id , Src: $src , src_num_INDEX: $src_num , src_num_LABEL: $src_label");
if ( $src_label =~ /PROPOSAL/i ){
info("DEBUG - LCMerge -......src_num_LABEL: PROPOSAL => Searching for an EPIC identification.");
my $data = readFITSTable(
file => $srcList
, extension => 'SRCLIST'
, colname => [qw( INDEX LABEL RA DEC RATE DELTA_DISP DELTA_XDSP PROCESS FLAG )]
);
my $delta_disp_proposal = $data->{DELTA_DISP}[0];
my $delta_xdsp_proposal = $data->{DELTA_XDSP}[0];
my %label_rate;
for my $i ( 2 .. $data->{-nrows} -1 ){
my $delta_disp = $data->{DELTA_DISP}[$i];
my $delta_xdsp = $data->{DELTA_XDSP}[$i];
my $label = $data->{LABEL}[$i];
my $process = $data->{PROCESS}[$i];
my $dist = 60 * (($delta_disp - $delta_disp_proposal)**2 + ($delta_xdsp - $delta_xdsp_proposal)**2)**(1/2);
my $dist_limit = 6;
if ($dist < $dist_limit and !$process){
$label_rate{$label} = $data->{RATE}[$i];
}
}
$epic_src_num = List::Util::reduce { $label_rate{$b} > $label_rate{$a} ? $b : $a } keys %label_rate;
$epic_src_num =~ s/EPIC0+//;
info("DEBUG - LCMerge -......EPIC identification, EPIC SRC ID = $epic_src_num");
} else {
$epic_src_num = $src_label;
$epic_src_num =~ s/EPIC//;
$epic_src_num = $epic_src_num + 0;
info("DEBUG - LCMerge -......EPIC identification, EPIC SRC ID = $epic_src_num");
}
if ( $inst =~ /rgs1/){
$rgs1{$src_num} = $epic_src_num;
} elsif ( $inst =~ /rgs2/){
$rgs2{$src_num} = $epic_src_num;
}
}
}
}
my @epic_sources;
my %rrgs1 = reverse %rgs1;
my %rrgs2 = reverse %rgs2;
info("******************************************************************************************");
info("DEBUG - LCMerge - RGS products vs EPIC identifications:");
foreach my $index ( sort ( keys %rgs1 ) ){
info("DEBUG - LCMerge -......rgs1 index: $index => EPIC SRC ID: $rgs1{$index}");
push @epic_sources, $rgs1{$index};}
info("DEBUG - LCMerge -......");
foreach my $index ( sort ( keys %rgs2 ) ){
info("DEBUG - LCMerge -......rgs2 index: $index => EPIC SRC ID: $rgs1{$index}");
push @epic_sources, $rgs1{$index};}
info("******************************************************************************************");
@epic_sources = grep { $_ ne '' } @epic_sources;
my @epic_sources_uniq = uniq @epic_sources;
my @rgs_keys = uniq (keys %rgs1, keys %rgs2);
info("******************************************************************************************");
info("********** DEBUG - LCMerge - Loop over RGS1/2 source products. RGS indices: @rgs_keys **********");
info("******************************************************************************************");
my $check_pn_timing = 0;
foreach my $rgs_key ( sort @rgs_keys ){
info("******************************************************************************************");
info("DEBUG - LCMerge -......RGS1/2 index: $rgs_key . EPIC SRC IDs: $rgs1{$rgs_key} , $rgs2{$rgs_key}");
info("******************************************************************************************");
my $check_rgs1_proposal = 0;
my $check_rgs2_proposal = 0;
if ( $rgs_key == 1 ){
$check_rgs1_proposal = 1;
$check_rgs2_proposal = 1;
}
my $lcurveFileList;
if ( $rgs_key == 1 && ! $rgs1{$rgs_key}){
info("DEBUG - LCMerge -......PROPOSAL + NO EPIC SRC ID");
info("DEBUG - LCMerge -......All RGS products will be merged only with PN TIMING TS.");
my %instr_timing;
my @lines = ();
my @epicTSList;
foreach my $instr (qw(epn emos1 emos2)){
@epicTSList = findFile(
class => 'product'
, instrument => $instr
, src_num => '0'
, content => 'EPIC source timeseries'
, format => 'FITS'
);
next unless ( scalar(@epicTSList) != 0 );
$instr_timing{$instr} = \@epicTSList;
info("DEBUG - LCMerge -......Instrument: $instr , Timing TS files: @{$instr_timing{$instr}}");
@{$instr_timing{$instr}} = sortByTSTART(@{$instr_timing{$instr}});
foreach my $i ( @{$instr_timing{$instr}} ){
chomp $i;
info("DEBUG - LCMerge -......$instr TSs: $i");
push (@lines, "$i\n");
}
}
if ( scalar(@lines) ){
$check_pn_timing = 1;
info("DEBUG - LCMerge -......RGS PROPOSAL products + EPIC Timing products.");
}
if ($check_rgs1_proposal || $check_rgs2_proposal){
my @rgs1TSList = findFile(
class => 'product'
, instrument => 'rgs1'
, src_num => '1'
, content => 'RGS SOURCE TIMESERIES'
, format => 'FITS'
);
info("DEBUG - LCMerge -......RGS PROPOSAL product RGS1 : @rgs1TSList");
my @rgs2TSList = findFile(
class => 'product'
, instrument => 'rgs2'
, src_num => '1'
, content => 'RGS SOURCE TIMESERIES'
, format => 'FITS'
);
info("DEBUG - LCMerge -......RGS PROPOSAL product RGS2 : @rgs2TSList");
push @rgsTSList, @rgs1TSList;
push @rgsTSList, @rgs2TSList;
@rgsTSList = grep { $_ ne '' } @rgsTSList;
@lines = grep { $_ ne '' } @lines;
if (@lines and scalar (grep { $_ !~ /\/\/\// } $lines[-1]) and @rgsTSList ) { push (@lines, "///\n"); }
@rgsTSList = sortByTSTART(@rgsTSList);
foreach my $i ( @rgsTSList ){
chomp $i;
info("DEBUG - LCMerge -......R1/2 LCs: $i");
push (@lines, "$i\n");
}
@lines = grep { $_ ne '' } @lines;
$src_files{1} = '0';
$src_files{2} = \@lines;
my $lcurvePlot_check = lcurvePlot(%src_files);
if ( $lcurvePlot_check ){ info("DEBUG - LCMerge - lcurvePlot completed for RGS PROPOSAL + NO EPIC identification");
info("******************************************************************************************");
info("DEBUG - LCMerge -"); }
}
} else {
info("DEBUG - LCMerge - RGS products with EPIC identification, EPIC SRC ID: @epic_sources_uniq");
foreach my $epic_src_num ( sort @epic_sources_uniq ){
my @pnTimingTSList;
if ( ! $check_pn_timing) {
@pnTimingTSList = findFile(
class => 'product'
, instrument => 'epn'
, src_num => '0'
, content => 'EPIC source timeseries'
, format => 'FITS'
);
} else {
info("DEBUG - LCMerge -......PN TIMING TS was included with RGS PROPOSAL.");
}
if ( scalar(@pnTimingTSList) ){
@epnTSList = @pnTimingTSList;
info("DEBUG - LCMerge -......PN TIMING TS: @epnTSList");
} else {
@epnTSList = findFile(
class => 'product'
, instrument => 'epn'
, src_num => $epic_src_num
, content => 'EPIC source timeseries'
, format => 'FITS'
);
info("DEBUG - LCMerge -......PN TS for EPIC SRC ID $epic_src_num: @epnTSList");
}
@emos1TSList = findFile(
class => 'product'
, instrument => 'emos1'
, src_num => $epic_src_num
, content => 'EPIC source timeseries'
, format => 'FITS'
);
@emos2TSList = findFile(
class => 'product'
, instrument => 'emos2'
, src_num => $epic_src_num
, content => 'EPIC source timeseries'
, format => 'FITS'
);
my $src_num_1 = $rrgs1{$epic_src_num};
if ( defined($src_num_1) ){
@rgs1TSList = findFile(
class => 'product'
, instrument => 'rgs1'
, src_num => $src_num_1
, content => 'RGS SOURCE TIMESERIES'
, format => 'FITS'
);
}
my $src_num_2 = $rrgs2{$epic_src_num};
if ( defined($src_num_2) ){
@rgs2TSList = findFile(
class => 'product'
, instrument => 'rgs2'
, src_num => $src_num_2
, content => 'RGS SOURCE TIMESERIES'
, format => 'FITS'
);
}
my @lines = ();
@epnTSList = sortByTSTART(@epnTSList);
foreach my $i ( @epnTSList ){
chomp $i;
info("DEBUG - LCMerge -......PN LCs for EPIC SRC ID $epic_src_num : $i");
push (@lines, "$i\n");
}
@lines = grep { $_ ne '' } @lines;
if (@lines and scalar (grep { $_ !~ /\/\/\// } $lines[-1]) and @emos1TSList ) { push (@lines, "///\n"); }
@emos1TSList = sortByTSTART(@emos1TSList);
foreach my $i ( @emos1TSList ){
chomp $i;
info("DEBUG - LCMerge -......M1 LCs for EPIC SRC ID $epic_src_num : $i");
push (@lines, "$i\n");
}
if (@lines and scalar (grep { $_ !~ /\/\/\// } $lines[-1]) and @emos2TSList ) { push (@lines, "///\n"); }
@emos2TSList = sortByTSTART(@emos2TSList);
foreach my $i ( @emos2TSList ){
chomp $i;
info("DEBUG - LCMerge -......M2 LCs for EPIC SRC ID $epic_src_num : $i");
push (@lines, "$i\n");
}
undef @rgsTSList;
push @rgsTSList, @rgs1TSList;
push @rgsTSList, @rgs2TSList;
@rgsTSList = grep { $_ ne '' } @rgsTSList;
if (@lines and scalar (grep { $_ !~ /\/\/\// } $lines[-1]) and @rgsTSList ) { push (@lines, "///\n"); }
@rgsTSList = sortByTSTART(@rgsTSList);
foreach my $i ( @rgsTSList ){
chomp $i;
info("DEBUG - LCMerge -......RGS1/2 LCs for RGS index: $src_num_1,$src_num_2 EPIC SRC ID $epic_src_num: $i");
push (@lines, "$i\n");
}
$src_files{1} = $epic_src_num;
$src_files{2} = \@lines;
my $lcurvePlot_check = lcurvePlot(%src_files);
if ( $lcurvePlot_check ){ info("DEBUG - LCMerge - lcurvePlot completed for RGS + EPIC SRC IDs: $epic_src_num");
info("******************************************************************************************");
info("DEBUG - LCMerge -"); }
}
}
}
info("******************************************************************************************");
info("********** DEBUG - LCMerge - Rest of TS's only for EPIC: PN + MOS's **********************");
info("******************************************************************************************");
my $summaryList = findFile(
class => "product"
, instrument => 'epic'
, content => 'EPIC summary source list'
, format => 'FITS'
);
if ( fileExists( file => $summaryList )){
my $filteredSummaryList = newFile(
class => 'intermediate'
, instrument => 'epic'
, content => 'Filtered EPIC summary source list'
);
doCommand(
'fselect'
, infile => $summaryList."[SRCLIST]"
, outfile => $filteredSummaryList
, expr => "TSERIES"
)
or return exception();
my $src_num_withTSeries = readFITSColumn(
file => $filteredSummaryList
, extension => "SRCLIST"
, column => "SRC_NUM"
);
for ( my $i=0 ; defined $$src_num_withTSeries[$i] ; $i++ ){
my $epic_src_num = $$src_num_withTSeries[$i];
if ( $epic_src_num ~~ @epic_sources_uniq ){ next; }
info("DEBUG - LCMerge - EPIC source with TS but not RGS product. EPIC SRC ID: $epic_src_num");
undef @epnTSList;
undef @emos1TSList;
undef @emos2TSList;
@epnTSList = findFile(
class => 'product'
, instrument => 'epn'
, src_num => $epic_src_num
, content => 'EPIC source timeseries'
, format => 'FITS'
);
@emos1TSList = findFile(
class => 'product'
, instrument => 'emos1'
, src_num => $epic_src_num
, content => 'EPIC source timeseries'
, format => 'FITS'
);
@emos2TSList = findFile(
class => 'product'
, instrument => 'emos2'
, src_num => $epic_src_num
, content => 'EPIC source timeseries'
, format => 'FITS'
);
my @lines = ();
@epnTSList = sortByTSTART(@epnTSList);
@emos1TSList = sortByTSTART(@emos1TSList);
@emos2TSList = sortByTSTART(@emos2TSList);
foreach my $i ( @epnTSList ){
chomp $i;
info("DEBUG - LCMerge -......PN LCs for EPIC SRC ID $epic_src_num : $i");
push (@lines, "$i\n");
}
@lines = grep { $_ ne '' } @lines;
if (@lines and scalar (grep { $_ !~ /\/\/\// } $lines[-1]) and @emos1TSList ) { push (@lines, "///\n"); }
foreach my $i ( @emos1TSList ){
chomp $i;
info("DEBUG - LCMerge -......M1 LCs for EPIC SRC ID $epic_src_num : $i");
push (@lines, "$i\n");
}
if (@lines and scalar (grep { $_ !~ /\/\/\// } $lines[-1]) and @emos2TSList ) { push (@lines, "///\n"); }
foreach my $i ( sort @emos2TSList ){
chomp $i;
info("DEBUG - LCMerge -......M2 LCs for EPIC SRC ID $epic_src_num : $i");
push (@lines, "$i\n");
}
$src_files{1} = $epic_src_num;
$src_files{2} = \@lines;
my $lcurvePlot_check = lcurvePlot(%src_files);
if ( $lcurvePlot_check ){ info("DEBUG - LCMerge - lcurvePlot completed for EPIC SRC ID: $epic_src_num");
info("******************************************************************************************");
info("DEBUG - LCMerge -"); }
}
}
return success();
}
sub sortByTSTART {
my @files = @_;
my (%file_tstart, @lines);
foreach my $i ( sort @files ){
my $tstart = readFITSKeyword(file => $i, extension => 'RATE', keyword => 'TSTART');
$file_tstart{$i} = $tstart;}
foreach my $i (sort { $file_tstart{$a} <=> $file_tstart{$b} } keys %file_tstart){
push (@lines, "$i\n");}
return @lines;
}
sub lcurvePlot {
my (%src_files) = @_;
my $src_num = $src_files{1};
my @lines = @{$src_files{2}};
my $inputSeries = 1 + scalar (grep { $_ =~ /\/\/\// } @lines);
info("******************************************************************************************");
info("DEBUG - LCMerge - Merged TS plot. Number of instruments/plots: $inputSeries");
my $lcurveFileList = newFile(
class => 'intermediate'
, src_num => $src_num
, content => 'lcurve fileList'
, format => 'ASCII'
);
if ( fileExists( file => $lcurveFileList )){unlink ( $lcurveFileList );}
writeASCIIFile( name => $lcurveFileList
, text => \@lines
, append => 1
);
my $lcurvePlot = newFile(
class => 'intermediate'
, src_num => $src_num
, content => 'lcurve plot'
, format => 'PS'
);
my $lcurveCommands = newFile(
class => 'intermediate'
, src_num => $src_num
, content => 'lcurve commands'
, format => 'ASCII'
);
if ( fileExists( file => $lcurveCommands )){unlink ( $lcurveCommands );}
my @commands = ();
my $j = 1;
my @files = grep { $_ !~ /\/\/\// } @lines;
my $files_size = scalar (@files);
my $number_of_RGS_LCs = scalar ( grep { $_ =~ /R1/ } @files) + scalar ( grep { $_ =~ /R2/ } @files);
if ( $files_size eq $number_of_RGS_LCs){
info("DEBUG - LCMerge - List RGS TSs to merge. No plot for those cases.");
return;
}
foreach my $file (@files){
info("DEBUG - LCMerge -......\@files : $file");
}
foreach my $i (@files){
my $label;
$j++;
if (scalar ( grep { $_ =~ /PN/ } $i)) {$label = "PN";}
if (scalar ( grep { $_ =~ /M1/ } $i)) {$label = "M1";}
if (scalar ( grep { $_ =~ /M2/ } $i)) {$label = "M2";}
if (scalar ( grep { $_ =~ /R1/ } $i)) {$label = "RGS";}
if (scalar ( grep { $_ =~ /R2/ } $i)) {$label = "RGS";}
info("DEBUG - LCMerge -......Plot label: $label");
push (@commands, "win $j\n");
push (@commands, "lab y $label\n");
if (scalar (grep { $_ =~ /PN/ } @commands) > 1) {
splice @commands, -2;
$j = $j -1;
}
if (scalar (grep { $_ =~ /M1/ } @commands) > 1) {
splice @commands, -2;
$j = $j -1;
}
if (scalar (grep { $_ =~ /M2/ } @commands) > 1) {
splice @commands, -2;
$j = $j -1;
}
if (scalar (grep { $_ =~ /RGS/ } @commands) > 1) {
splice @commands, -2;
$j = $j -1;
}
}
push (@commands, "Device $lcurvePlot/vps\n");
push (@commands, "plot\n");
push (@commands, "Device\n");
push (@commands, "exit\n");
writeASCIIFile( name => $lcurveCommands
, text => \@commands
, append => 1
);
my $command = 'lcurve nser=' . $inputSeries . ' cfile1="@' . $lcurveFileList .
'" window="-" plot=yes plotdev="/null" dtnb="INDEF" nbint="INDEF" outfile=" " plotdnum=' .
$inputSeries .' < ' . $lcurveCommands;
doCommand("$command");
if ( fileExists(file => $lcurvePlot)){
info("DEBUG - LCMerge - Merged timeseries created");
my $lcurveInterPDFPlot = newFile(
class => 'intermediate'
, instrument => 'epic'
, src_num => $src_num
, content => 'EPIC source lcurve timeseries plot'
, format => 'PDF'
);
PStoPDF(
source => $lcurvePlot
, destination => $lcurveInterPDFPlot
) or return exception();
} else {
info("DEBUG - LCMerge - Warning: Merged timeseries was not created");
}
my $LCMergePDFplot = newFile(
class => 'product'
, instrument => 'epic'
, src_num => $src_num
, content => 'EPIC source timeseries plot'
, format => 'PDF'
);
doCommand(
'plotTimeSeries.py'
, infile => $lcurveFileList
, outfile => $LCMergePDFplot
, outformat => 'pdf'
);
}
1;