package EPICSources;
use strict;
use English;
use Carp;
use vars qw(@ISA $VERSION $name $author $date $version @instList $numberOfStreams);
@ISA = qw(Module);
use ModuleResources;
use ModuleUtil;
$name = __PACKAGE__;
$VERSION = '0.17';
$version = $VERSION;
$author = 'Anja Schroeder,Ian Stewart,Duncan John Fyfe, Duncan Law-Green, Jose V Perea';
$date = '2020-11-27';
#
# ChangeLog
# =========
#
# version 0.17 - 2020-11-27 (JVP)
# ------------
#
# Skip EPIC-OM SED plot if there is not OM src list or EPIC src list
# Skip also if there are no OM counterparts for EPIC sources
#
# version 0.16 - 2020-10-05 (JVP)
# ------------
#
# Include EPIC-OM combined Spectral Energy Distribution plot (SED)
# EPIC-OM sources cross-match, EPIC and OM finding charts and SED plot
# New Python task: epicSED.py
#
# version 0.15 - 2014-02-28 (EOP)
# ------------
#
# + Following the recommendation in the SAS Users Guide and SAS threads, MOS
# spectral products should be created after applying #XMMEA_EM flag
# filtering.
# NOTE: These flag masks are applied the same in EPICSourceProducts.pm
#
# version 0.14 - 2012-03-15 DLG
# ------------
#
# + esources mintotalcts 100.
#
# version 0.14 - 2012-03-15 DLG
# ------------
#
# + esources mintotalcts 300. Again.
#
# version 0.13 - 2012-03-13 DLG
# ------------
#
# + try reverting esources mintotalcts to 500 to avoid xspec spectral extraction failures
#
# version 0.12 - 2012-01-19 DLG
# ------------
#
# + Added Attitude GTI filtering to $specExpression and $timeExpression to
# propagate to SSP spectra and time series
#
# version 0.11 - 2012-01-19 DLG
# ------------
#
# + Altered mintotalcts in esources from 500 to 300 to test optimised
# faint source extraction
#
# version 0.10 - 2009-09-14 DLG
# ------------
#
# + Added test for MOS W6 and PN SW modes. Exclude these from SSP generation.
#
# version 0.09 - 2006-08-11 DJF
# ------------
#
# + Added a test for PN CCD in low gain mode. Exclude PN exposures where some CCDs are in low gain mode.
#
# version 0.08 - 2006-03-03 ACS
# ------------
# + A couple of corrections to the perl code introduced by the changes the day before.
#
# version 0.07 - 2006-03-02 ACS
# ------------
# + Added determination of minimum TSTART and maximum TSTOP of all valid exposures
# and inserted them into each filtered event list for light curves
# + Commented out selection by data mode ..Prime.. since FU imaging event lists
# did not get filtered event lists.
#
# version 0.06 - 2005-12-02 IMS
# ------------
# + Now really added --setflags=yes (d'oh!!)
#
# version 0.05 - 2005-12-01 IMS
# ------------
# + Added --setflags=yes to esources call. This sets the SPECTRA and TSERIES flag columns of the srcmatch list to T for those sources selected for SSPs.
#
# version 0.04 - 2005-11-29 IMS
# ------------
# + Changed the flag mask selections to those agreed upon.
#
# version 0.03 - 2005-10-14 IMS
# ------------
# + Shortened the content string for the SSP source list
#
# version 0.02 - 2005-09-08 RGW
# ------------
#
# + exit gracefully if required source list doesn't exist
#
# version 0.01 - 2005-06-09 IMS
# ------------
#
# + First draft - ported from ACS script.
#
# Declare instruments of interest
@instList = qw(all);
# Number of Streams
sub numberOfStreams
{
return 1;
}
# Rules evaluated to determine when this modules can run.
sub evaluateRules
{
# If upstream module has been ignored, so if this one
return ignore()
if ignored(
module => 'SrcMerge'
, instrument => 'all'
, stream => 1
);
#
return start()
if complete(
module => 'SrcMerge'
, instrument => 'all'
, stream => 1
)
}
# The guts of the module
sub performAction
{
info("Module version number: $version");
# First fetch list of all suitable EPIC event lists (imaging, not cal. closed)
#
my (@eventListSets, $content, $timeExpression, $specExpression);
my (@filteventListSets, $tlmin, $tlmax, @tmin, @tmax, @tmpmin, @tmpmax);
# get attitude GTI file
my $attGTIFile = findFile(
class => 'intermediate'
, instrument => 'all'
, content => 'Attitude GTI'
);
# add attitde GTI filtering if file exists
my $attGTIExpr = "";
$attGTIExpr = "&>I($attGTIFile,TIME)" if $attGTIFile;
foreach my $inst ( "emos1", "emos2", "epn" )
{
if ($inst eq 'epn')
{
$content= "EPIC PN imaging mode event list";
$timeExpression
# ="((PATTERN<=4)&&((FLAG & 0x2fb002c)==0)"
="((PATTERN<=4)&&((FLAG & 0xfffffef)==0)"
."&&(PI>=200)&&(PI<=12000))"
.$attGTIExpr;
$specExpression
="((PATTERN<=4)&&(FLAG==0)&&(PI>=200)&&(PI<=12000))"
.$attGTIExpr;
}
else # mos
{
$content= "EPIC MOS imaging mode event list";
$timeExpression
# ="((PATTERN<=12)&&((FLAG & 0x766b0000)==0)"
="((PATTERN<=12)&&((FLAG & 0x766ba000)==0)"
."&&(PI>=200)&&(PI<=12000))"
.$attGTIExpr;
# $specExpression = $timeExpression;
$specExpression
# ="((PATTERN<=12)&&((FLAG & 0xfffffeff)==0)"
="((PATTERN<=12)&&((FLAG & 0x766ba000)==0)"
."&&(PI>=200)&&(PI<=12000))"
.$attGTIExpr;
}
foreach my $exp_id ( listExposures( instrument => $inst ) )
{
# Find imaging event file for this exposure
my $evFile=findFile(
class => 'product'
, instrument => $inst
, exp_id => $exp_id
, content => $content
);
next if (! $evFile );
if ($inst eq 'epn')
## Select EPIC PN
{
my $lowgain = {};
## Exclude PN CCD low gain mode from SSP generation
if (my @lowgain = ModuleUtil::pnLowGainWarn($evFile , $lowgain) )
{
info( 'Some CCD in low gain mode :' , join (' ,' , @lowgain) );
info( 'Excluding this exposure from source detection.' );
next;
}
## Exclude PN SW mode from SSP generation
my $submode = readFITSKeyword(
file => $evFile
, extension => 'EVENTS'
, keyword => 'SUBMODE'
);
# EC (14 Nov 2013): don't skip, at request of PR
#if ($submode =~ /PrimeSmallWindow/) {
# info ( 'PN PrimeSmallWindow mode detected: excluding from SSP. ');
# next;
#}
} else {
## Select EPIC MOS
## Exclude MOS W6 mode from SSP generation
my $submode = readFITSKeyword(
file => $evFile
, extension => 'EVENTS'
, keyword => 'SUBMODE'
);
if ($submode =~ /PrimePartialW6/) {
info ( 'MOS PrimePartialW6 mode detected: excluding from SSP. ');
next;
}
}
# Check instrument filter
my $filter = getExposureProperty(
instrument => $inst
, exp_id => $exp_id
, name => 'filter'
);
next if $filter =~ /Cal/;
# Derive Tstart and Tstop times for time series
$tlmin = readFITSKeyword(
file => $evFile
, extension => 'EVENTS'
, keyword => 'TSTART'
);
$tlmax = readFITSKeyword(
file => $evFile
, extension => 'EVENTS'
, keyword => 'TSTOP'
);
# Add to lists
push @eventListSets, $evFile;
push @tmin, $tlmin;
push @tmax, $tlmax;
my $forTimeFilteredEvList = newFile(
class => 'intermediate'
, instrument => $inst
, exp_id => $exp_id
, content => 'Event list for light curves'
);
doCommand(
'evselect'
, table => "${evFile}:EVENTS"
, keepfilteroutput => 'yes'
, withfilteredset => 'yes'
, filteredset => $forTimeFilteredEvList
, writedss => 'yes'
, updateexposure => 'yes'
, expression => $timeExpression
, destruct => 'yes'
)
or return exception();
push @filteventListSets, $forTimeFilteredEvList;
my $forSpecFilteredEvList = newFile(
class => 'intermediate'
, instrument => $inst
, exp_id => $exp_id
, content => 'Event list for spectra'
);
doCommand(
'evselect'
, table => "${evFile}:EVENTS"
, keepfilteroutput => 'yes'
, withfilteredset => 'yes'
, filteredset => $forSpecFilteredEvList
, writedss => 'yes'
, updateexposure => 'yes'
, expression => $specExpression
, destruct => 'yes'
)
or return exception();
} # end loop over $exp_id
} # end loop over $inst
# Return immediately if there are no suitable exposures.
unless (@eventListSets)
{
info("No suitable exposures found");
return success();
}
# Determine min and max of all tstart and tstop times
@tmpmin = sort(@tmin);
@tmpmax = sort(@tmax);
my $tminall = $tmpmin[0];
my $tmaxall = $tmpmax[$#tmpmax];
# Add TLMIN1 and TLMAX1 keywords
doCommand('addattribute'
,set => [ @filteventListSets ]
,attributename => [qw(EVENTS:TLMIN1 EVENTS:TLMAX1)]
,attributetype => ['real ' , 'real']
,realvalue => [ $tminall , $tmaxall ]
,attributecomment => [ 'Time series start' , 'Time series stop' ]
);
# Generate Source list:
my $srcmatchSrcFile = findFile(
class => 'product'
, instrument => 'epic'
, content => 'EPIC summary source list'
, format => 'FITS'
);
return success() unless $srcmatchSrcFile;
my $sspSrcFile = newFile(
class => 'intermediate'
, instrument => 'epic'
# , content => 'Source list for source-specific products'
, content => 'SSP source list'
);
doCommand(
'esources'
, eventsets => [@eventListSets]
, srctab => "${srcmatchSrcFile}:SRCLIST"
, outsrcset => $sspSrcFile
, mintotalcts => 100
, mindetml => 15
, minmaskfrac => 0.5
, numhexchar => 3
, setflags => 'yes'
)
or return exception();
#########################################################################################
# EPIC - OM combined Spectral Energy Distribution (SED)
info("Start producing EPIC-OM SED plot");
# Find OM and EPIC source lists
# Get om summary source list, if it exists, and has non-zero rows <<==== PENDING
my $OMsrcList = findFile(
class => 'product'
, instrument => 'om'
, content => 'OM OBSERVATION SOURCE LIST'
);
my $EPICsrcList = findFile(
class => "product"
, instrument => 'epic'
, content => 'EPIC summary source list'
, format => 'FITS'
);
info("SED DEBUG -- EPICsrcList : $EPICsrcList");
info("SED DEBUG -- OMsrcList : $OMsrcList");
if ( !fileExists( file => $EPICsrcList) || !fileExists( file => $OMsrcList ) ){
info("EPICsrcList or OMsrcList not found. Skipping EPIC-OM SED plot production.");
return success();
}
# Get EPIC sources with OM counterpart:
# SRC_NUM entries with column N_OPT_CO != -1
#
my @sources = nonNullCol( $EPICsrcList, 'SRCLIST' );
info("SED DEBUG -- EPIC sources with OM counterpart: @sources");
if ( ! @sources ){
info("No OM counterparts for EPIC sources. Skipping EPIC-OM SED plot production.");
return success();
}
my (@sourcesFiles, $srcFile, $outSEDplot);
foreach my $src_num (@sources){
$outSEDplot = newFile( class => 'product'
, instrument => 'epic'
, src_num => $src_num
, content => "EPIC-OM SOURCE SED PLOT"
, format => 'PNG'
);
$srcFile = "$src_num:$outSEDplot";
push (@sourcesFiles, $srcFile);
}
my $sourcesFiles = join(',', @sourcesFiles);
info("sourcesFiles : $sourcesFiles");
# EPIC 3Colour image
my $epic3COLimage = findFile(
class => 'product'
, instrument => 'epic'
, band => 8
, content => 'EPIC THREECOLOUR IMAGE'
, format => 'FITS'
);
info("SED DEBUG -- epic3COLimage : $epic3COLimage");
# OM images from modes:
my @modes = (
"OM USER WINDOWS SKY IMAGE MOSAIC"
, "OM RUDI-5 SKY IMAGE MOSAIC"
, "OM FULL-FRAME LORES SKY IMAGE MOSAIC"
, "OM FULL-FRAME HIRES SKY IMAGE MOSAIC"
);
my @omimages;
foreach my $mode (@modes){
my @images = findFile(
class => 'product'
, instrument => 'om'
, content => "$mode"
, format => 'FITS'
);
push (@omimages, @images);
}
my $omimages = join(',', @omimages);
info("SED DEBUG -- omimages : @omimages");
#
print "Start epicSED.py ...\n";
doCommand(
'epicSED.py'
, omsrclist => $OMsrcList
, epicsrclist => $EPICsrcList
, epic3COLimage => $epic3COLimage
, omimages => $omimages
, plotformat => 'png'
, outfiles => $sourcesFiles
)
or return exception();
return success();
}
# Sub for testing. To be moved to ModuleUtils.pm
sub nonNullCol
{
my ($file, $ext, $colname, $cols) = @_;
my $data = readFITSTable(
file => $file
, extension => $ext
, colname => [qw( SRC_NUM N_OPT_CO )]
#, colname => [$cols]
);
my @srcs;
for my $i ( 0 .. $data->{-nrows} -1 ){
my $src_num = $data->{SRC_NUM}[$i];
my $n_opt_co = $data->{N_OPT_CO}[$i];
push(@srcs, $src_num) unless ($n_opt_co < 0);
}
return @srcs;
}
1;