package RGSProducts;
use strict;
use English;
use Carp;
use vars
qw(@ISA $VERSION $name $author $date $version @instList $numberOfStreams);
@ISA = qw(Module);
use ModuleResources;
use ModuleUtil;
use Astro::FITS::CFITSIO qw(TRUE FALSE);
use Regexp::Common qw(boolean);
use POSIX; # for floor
# Declare identity, version, author, date, etc.
$name = __PACKAGE__;
$VERSION = '2.34';
$version = $VERSION;
$author = 'Richard West,Dean Hinshaw,Duncan John Fyfe,Ian Stewart,Duncan Law-Green,Eduardo Ojero,Ed Chapin, Jose V Perea, Laura Tomas';
$date = '2017-12-14';
# ChangeLog
# =========
#
# version 2.34 - 2017-12-14 LTC
# ------------
#
# + RGS response matrices for order 2 spectra
#
# version 2.33 - 2017-11-29 JVP
# ------------
#
# + Skip source Time Series plot ('elcplot') if the input FITS file has not BACKV column (bkg rate) or
# if that column has some negative value. (All BACKV values have to be > 0)
#
# version 2.32 - 2016-12-16 JVP
# ------------
#
# + No spectral binning in RGS spectra ( rgsspecplot rebin = Y -> N )
#
# version 2.31 - 2016-09-19 JVP
# ------------
#
# + Time binsize for RGS LCs calculated from the total counts in the spectrum and the total exposure time.
# Keeping at least 20 counts per bin.
#
# version 2.30 - 2016-08-02 JVP
# ------------
#
# + No yrange limits in RGS Image if observation is processed in the target reference frame ( SSO: --sso_object )
# No limits in the 'Cross dispersion angle' plot image.
#
# version 2.29 - 2016-03-29 LT
# ------------
#
# + New final product: RGS source timeseries plot (pdf)
#
# version 2.28 - 2014-10-01 JVP
# ------------
#
# + rgsbkgmodelTest corrected (rgsbkgmodel-1.4.1) to exit without exception (error code = 0 now).
# Command output captured from stdout to test if 'rgsbkgmodel' can be passed.
#
# version 2.27 - 2014-02-13 EC
# ------------
#
# + No rgsbkgmodel if SpectroscopySmallWindow mode (Mantis #0007175)
#
# version 2.26 - 2013-10-02 EC
# ------------
#
# + Switch to LAMBDA from BETA_CORR, including updates to ximagemin, and ximagemax (Mantis 0007141)
#
# Version 2.25 - 2013-06-11 EC
# ------------
#
# + Merge changes from SOC version.
#
# version 2.24.1 - 2011-08-09 EOP
# --------------
#
# + Creation of GIF files by implot replaced by direct PNG file writing.
# This same code is used on both 32 and 64-bit pipelines.
#
# version 2.24 2011-02-10 DLG
# ------------
#
# + Disable b/g subtraction for light curve generation in SpectroscopySmallWindow
#
# version 2.23 2011-02-07 DLG
# ------------
#
# + Enable wavelength binning for whole-field spectrum
#
# version 2.22 2011-02-07 DLG
# ------------
#
# + Re-enable RGS wavelength binning, timeseries generation
#
# version 2.21 2011-02-07 DLG
# ------------
#
# + Eliminate use of GIF bitmap images
#
# version 2.20 2011-02-07 DLG
# ------------
#
# + Disable background generation for SpectroscopySmallWindow mode
#
# version 2.19 2010-09-07 DLG
# ------------
#
# + Disable rgslccorr at SRR request
#
# version 2.18 2010-06-02 DLG
# ------------
#
# + Set value for rgsspectrum spectralbinning to 'beta' explicitly
#
# version 2.17 2009-10-06 DLG
# ------------
#
# + Changed $timebinsize from 500 to 1000 on AP's advice
#
# version 2.16 2009-10-05 DLG
# ------------
#
# + Revert changes to spectral wavelength binning
# await fix to SAS task rgsspecplot
#
# version 2.15 2009-09-24 DLG
# ------------
#
# + Test implementation of spectral wavelength binning
#
# version 2.14 2009-05-20 DLG
# ------------
#
# + Added sourceid parameter to rgslccorr call
#
# version 2.13 2009-05-19 DLG
# ------------
#
# + Added test call to rgslccorr to evaluate RGS light curve generation in pipeline
#
# version 2.12 2007-01-23 DJF
# ------------
#
# + Read RGS source list for 'processable' source parameters and add these to the output FIT files.
# + Use common RGS source selection routine to select sources from which spectra will be created.
# + Write EPIC source parameters to RGS spectra and background when appropriate.
# + Fixed identification of the wholefield spectrum. FITS file row number and RGS source INDEX number were being used interchangeably. This is not always true.
#
# version 2.11 2006-09-20 DJF
# ------------
#
# + Record use of exposures for source specific product creation in the database.
#
# version 2.10 - 2006-07-01 DJF
# ------------
#
# + Updated to take account of parameter changes to rgsproc-1.24.1
#
# version 2.09 - 2006-01-19 DJF
# ------------
#
# + Added new sas task rgsbkgmodelTest to determine if rgsbkgmodel should be run
#
# version 2.08 - 2005-11-24 IMS
# ------------
# + Shortened name of whole-field intermediate ps file.
# + Whole-field spectrum should not be background-subtracted. Found that after the changes of the preceding version, I was making 2 products (with different names and content strings) which fulfilled this criterion! I've now removed the first call, and sent the fits output of the second to rgsspecplot.
#
# version 2.07 - 2005-11-15 IMS
# ------------
# + Modified call to rgsspectrum which creates the whole-field spectrum in response to Andy Pollock's suggestions on 8 Nov.
#
# version 2.06 - 2005-11-14 IMS
# ------------
#
# + Removed call to rgsfluxer - this is now done in module RGSMakeFluxed. Flag $CreateRgsFluxedSpectrum has also been disabled.
# + Initialized $i to 0 in for ($i etc) calls.
#
# version 2.05 - 2005-11-12 DJF
# ------------
#
# + rgsfluxer wasn't raising an exception on failure.
#
# version 2.04 - 2005-10-22 RGS
# ------------
#
# + fixed bug incorrectly numbering the source for the whole field spectrum
#
# version 2.03 - 2005-10-17 IMS
# ------------
#
# + Added creation of whole-field spectral products.
# + RGS source list was being 'findFiled' twice. I removed the redundant occurrence.
#
# version 2.02 - 2005-08-16 RGW
# ------------
#
# + name change on parameter to rgsbkgmodel
#
# version 2.01 - 2005-08-08 IMS
# ------------
#
# + Added call to rgsbkgmodel (after advice from A Pollock 5 Aug 05).
# + Response matrix and fluxed spectrum are now being made (as new products) (after advice from A Pollock 5 Aug 05).
### NOTE however that at this point in time final advice has not yet been received from SOC about how to make all the new products they want.
#
# version 2.00 - 2005-05-09 DJF
# ------------
#
# + Add explicit perl module headers. Previously these were supplied
# at compile time. This will make debugging and extending the modules
# through additional perl libraries easier.
#
# + Code layout made more uniform with perltidy
# + Standerdized date format (ccyy-mm-dd)
# + Standerdized author list to use comma separator
# + Make use of perl $VERSION magic. Now $Version = version = ''
##
# Version 1.34 - 2005-05-09 (DJF)
# ------------
#
# + Check that we have some events left after filtering.
#
# Version 1.33 - 22-Jan-2004 (DJF)
# ------------
#
# + rgsrmfgen rmfile should have been rmfset
#
# Version 1.32 - 21-Jan-2004 (DJF)
# ------------
#
# + Changed rgsrmfgen parameter file to rmfile to match change from SAS 5.x to 6
#
# Version 1.31 - 10-Dec-2003 (DJF)
# ------------
#
# + Adapted doCommand to use anonymous lists for list parameters.
#
# Version 1.30 - 17-July-2003 (DJF)
# ------------
#
# + Replace missing } else { in conditional selection of rgsimplot parameters.
#
# Version 1.29 - 25-June-2003 (DJF)
# ------------
#
# + Turned off production of fluxes spectra as intermediate products
#
# Version 1.28 - 25-June-2003 (DJF)
# ------------
#
# + Turn rgs response matrix from intermediate to final product.
#
# Version 1.27 - 08-April-2003 (DJF)
# ------------
#
# + Make rgs response matrices
#
# Version 1.26 - 26-July-2002 (DJF)
# ------------
#
# + Clear spectra list after each source.
#
# Version 1.25 - 25-July-2002 (DJF)
# ------------
#
# + Turned rgs total source spectra from intermediate files to products.
#
# Version 1.24 - 22-July-2002 (DJF)
# ------------
#
# + RGS response matrix and fluxed spectra removed.
# The SAS tasks will not be sufficiently mature until the
# next major SAS release.
# + Fixes for total spectrum creation
#
# Version 1.23 - 16-July-2002 (DJF)
# ------------
#
# + Fixed intermediate filename for new total spectrum
#
# Version 1.22 - 09-July-2002 (DJF)
# ------------
#
# + Added production of total source spectra
#
# Version 1.21 - 12-June-2002 (DJF)
# ------------
#
# + order is not a valid parameter for intermediate files
#
# Version 1.20 - 12-June-2002 (DJF)
# ------------
#
# + Bug fix, rgsorder is not a valid parameter for intermediate files
#
# Version 1.19 - 11-June-2002 (DJF)
# ------------
#
# + Added sas tasks rgsrmfgen and rgsfluxer
#
# Version 1.18 - 19-Feb-2002 (DH)
# ------------
#
# + Add rebin and mincounts parameters to rgsspecplot.Set mincounts parameter
# to 20.
#
# Version 1.17 - 30-Nov-2001 (DH)
# ------------
#
# + Remove limit on the rgs source number that it has to be less than 100.
# + Bring in line with rgsproc 1.3, specifically change to new rgsspectrum
# interface.
#
# Version 1.16 - 30-Nov-2001 (DH)
# ------------
#
# + For newFile call for creating a spectrum, need to specify
# both the order and the rgsorder attributes.Do the same
# for the background spectrum.
#
# Version 1.15 - 29-Nov-2001 (DH)
# ------------
#
# + Change order attribute to rgsorder in the newFile call
# for creating a spectrum.This is what is hould have been
# all along.
#
# Version 1.14 - 06-Nov-2001 (DJF)
# -------------------------------
#
# + remove duplicate declerations of $gifFile and $pngFile
#
# Version 1.13 - 06-Nov-2001 (DJF)
# -------------------------------
#
# + Change orientation of pdf plot (device =/PS -> /VPS).
#
# Version 1.10 - 1-Oct-2001 (DH)
# ------------
#
# + Replace dslatts calls with the hasFITSExtension function.
#
# Version 1.09 - 2001-05-22 (DJF)
# ------------
#
# + Added updateexposure=no to PI filtering evselect
# to eliminate an unecessary warning message.
#
# Version 1.08 - 2001-05-15 (DJF)
# ------------
#
# + Moved return when no event list is found
#
# Version 1.07 - 2001-05-15 (DJF)
# ------------
#
# + Add PI filtering of RGS event list.
#
# Version 1.06 - 2001-05-14 (DH)
# ------------
#
# + Incoporate changes to RGS event list.SOURCES table is now
# called SRCLIST, and SRC_SELECT column now called PROCESS.
#
# Version 1.05 - 2001-03-20 (DH)
# ------------
#
# + Better checks for presence of region extensions.Do
# not produce some plots if regions are missing.
#
# Version 1.04 - 2001-03-16 (DH)
# ------------
#
# + Print statement in version 1.03 was in the wrong place
#
# Version 1.03 - 2001-03-16 (DH)
# ------------
#
# + Print out version number in performAction() for
# tracking purposes.
#
# Version 1.02 - 2001-03-15 (DH)
# ------------
#
# + Check for existence of region extensions before producing
# spectral products for a given source.
#
# Version 1.01 - 2000-12-11 (DH)
# ------------
#
# + First production version.
#
# Declare list of instruments this module is interested in
@instList = qw(rgs1 rgs2);
# Number of streams
sub numberOfStreams
{
return numberOfExposures();
}
# Rules method
sub evaluateRules
{
# If upstream module has been ignored, so if this one
return ignore()
if ignored( module => 'ExpChooser'
, instrument => thisInstrument
, stream => thisStream
)
or ignored( module => 'GlobalHK'
, instrument => 'all'
, stream => 1
)
or ignored( module => 'RGSEvents'
, instrument => thisInstrument
, stream => thisStream
);
#
return start()
if complete( module => 'RGSEvents'
, instrument => thisInstrument
, stream => thisStream
);
}
# Action method
sub performAction
{
info("Module version number: $version");
# Work out which exposure this stream maps on to
# I added Response matrices and fluxed spectra but they have been removed for the moment.
# Set these flags to 1 to re-enable.
# my ( $CreateRgsResponseMatrix, $CreateRgsFluxedSpectrum )
# = ( 1, 1 );
# $CreateRgsResponseMatrix = 1
# if ( $CreateRgsFluxedSpectrum && !$CreateRgsResponseMatrix );
my $CreateRgsResponseMatrix = 1;
my $exp_id = exposureID( instrument => thisInstrument
, stream => thisStream
);
info("Exposure ID: $exp_id");
# Find event file products
my $evFile0 = findFile( class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS EVENT LIST'
);
# Find observation mode
my $mode = getExposureProperty(
exp_id => $exp_id
, name => 'mode'
);
# Return immediately if it doesn't exist (for whatever reason)
# or it is empty
return success() unless $evFile0;
return success() unless (fileExists( file => $evFile0 ));
return success unless ( numberFITSRows( file => $evFile0
, extension => 'EVENTS'
) > 0 );
my $evFile = newFile( class => 'intermediate'
, instrument => thisInstrument
, content => 'PI Filtered Event List'
);
my $eventfilter = "PI>170";
doCommand( 'evselect'
, table => "$evFile0:EVENTS"
, withfilteredset => 'y'
, keepfilteroutput => 'y'
, updateexposure => 'no'
, filteredset => $evFile
, expression => "PI > 170"
, writedss => 'no'
, destruct => 'yes'
) or return exception();
# Create RGS banana plot
my $bananaFile = newFile( class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS ENERGY-DISPERSION IMAGE'
);
doCommand( 'evselect'
, table => "$evFile:EVENTS"
, imageset => $bananaFile
, withimageset => 'yes'
, xcolumn => 'M_LAMBDA'
, withxranges => 'yes'
, ximagemin => 3.5
, ximagemax => 40
, ycolumn => 'PI'
, withyranges => 'yes'
, yimagemin => 0
, yimagemax => 3500
, updateexposure => 'yes'
, writedss => 'yes'
, withimagedatatype => 'yes'
, imagedatatype => 'Int32'
) or return exception();
# Find source (if available)
my @srcListPars;
my $srcList = findFile( class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS source list'
);
my %src_param;
my @selectedSources;
my $inst = uc thisInstrument;
if ($srcList)
{
%src_param = &ModuleUtil::get_rgssources_param( $srcList , 'SRCLIST' , thisInstrument );
@selectedSources = sort ( keys ( %src_param ) );
@srcListPars = (
withendispregionsets => 'yes'
, srclistset => $srcList
, srcidlist => [ @selectedSources ]
);
}
else
{
@srcListPars = ( withendispregionsets => 'no' );
}
# Create preview equivalent
my $gifFile = newFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS ENERGY-DISPERSION IMAGE'
, format => 'GIF'
);
my $pngFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS ENERGY-DISPERSION IMAGE'
, format => 'PNG'
);
doCommand(
'rgsimplot', endispset => $bananaFile
, withendispset => 'yes'
, withspatialregionsets => 'no'
, withspatialset => 'no'
, colour => 3
, device => '/PNG'
, plotfile => $pngFile
, orderlist => '1 2'
, @srcListPars
)
or return exception();
# my $pngFile = newFile(
# class => 'product'
# , instrument => thisInstrument
# , exp_id => $exp_id
# , content => 'RGS ENERGY-DISPERSION IMAGE'
# , format => 'PNG'
# );
# # Convert from GIF to PNG
# GIFtoPNG(
# source => $gifFile
# , destination => $pngFile
# )
# or return exception();
# Create camera image
my $cameraFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS IMAGE'
);
# No yrange limits if observation is processed in the target reference frame ( SSO: --sso_object )
#
my $sso_object = $ENV{PCMS_SSO_OBJECT};
if ( $sso_object ){
info("Observation is processed in the target reference frame ( SSO: --sso_object ). No yrange limits in RGS cross dispersion angle plot image.");
doCommand(
'evselect', table => "$evFile:EVENTS"
, imageset => $cameraFile
, withimageset => 'yes'
, xcolumn => 'M_LAMBDA'
, withxranges => 'yes'
, ximagemin => 3.5
, ximagemax => 40
, ycolumn => 'XDSP_CORR'
, withyranges => 'no' # <--
#, yimagemin => -7.e-4
#, yimagemax => 7.e-4
, updateexposure => 'yes'
, writedss => 'yes'
, withimagedatatype => 'yes'
, imagedatatype => 'Int32'
)
or return exception();
} else {
doCommand(
'evselect', table => "$evFile:EVENTS"
, imageset => $cameraFile
, withimageset => 'yes'
, xcolumn => 'M_LAMBDA'
, withxranges => 'yes'
, ximagemin => 3.5
, ximagemax => 40
, ycolumn => 'XDSP_CORR'
, withyranges => 'yes'
, yimagemin => -7.e-4
, yimagemax => 7.e-4
, updateexposure => 'yes'
, writedss => 'yes'
, withimagedatatype => 'yes'
, imagedatatype => 'Int32'
)
or return exception();
}
# Do not continue unless there is a background region
return success()
unless hasFITSExtension(
file => $srcList
, extension => "${inst}_BACKGROUND"
);
# Create preview equivalent
$gifFile = newFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS IMAGE'
, format => 'GIF'
);
$pngFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS IMAGE'
, format => 'PNG'
);
if ($srcList)
{
@srcListPars = (
withspatialregionsets => 'yes'
, srclistset => $srcList
, srcidlist => join( " ", @selectedSources )
);
}
else
{
@srcListPars = ( withspatialregionsets => 'no' );
}
doCommand(
'rgsimplot', withspatialset => 'yes'
, spatialset => $cameraFile
, withendispset => 'no'
, withendispregionsets => 'no'
, colour => 3
, device => '/PNG'
, plotfile => $pngFile
, orderlist => '1 2'
, @srcListPars
)
or return exception();
# # Convert from GIF to PNG
# GIFtoPNG(
# source => $gifFile
# , destination => $pngFile
# )
# or return exception();
my @alltotspectra;
my @matrices;
foreach my $srcnum (@selectedSources)
{
setExposureProperty( instrument => thisInstrument , exp_id => $exp_id , name => 'ssp' , value => 0+TRUE );
my $spregion = "$srcList:${inst}_SRC${srcnum}_SPATIAL";
my @spectra;
my @totspectra;
my $total_spec_counts_1st_order = 0.0;
my $exposure;
# Create background and source spectrum for each order
foreach my $order ( 1, 2 )
{
my $bkgSpecFile = '';
unless ($mode eq 'SpectroscopySmallWindow') {
$bkgSpecFile = newFile(
class => 'product'
, instrument => thisInstrument
, src_num => $srcnum
, exp_id => $exp_id
, order => $order
, rgsorder => $order
, content => 'RGS SOURCE BACKGROUND SPECTRUM'
);
};
my $srcSpecFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, order => $order
, rgsorder => $order
, content => 'RGS SOURCE SPECTRUM'
);
# Check for RGS Small window mode
# Skip background generation if found
if ($mode eq 'SpectroscopySmallWindow') {
info('SmallWindow mode detected. Skipping background generation...');
doCommand(
'rgsspectrum', evlist => $evFile
, srclist => $srcList
, source => $srcnum
, rebin => 1
, order => $order
, edgechannels => 2
, withspectrum => 'Y'
, spectrumset => $srcSpecFile
, bkgcorrect => 'N'
, withbkgset => 'N'
, spectrumbinning => 'lambda'
)
or return exception();
} else {
doCommand(
'rgsspectrum', evlist => $evFile
, srclist => $srcList
, source => $srcnum
, rebin => 1
, order => $order
, edgechannels => 2
, withspectrum => 'Y'
, spectrumset => $srcSpecFile
, bkgcorrect => 'Y'
, withbkgset => 'Y'
, spectrumbinning => 'lambda'
, bkgset => $bkgSpecFile
)
or return exception();
};
# Copy the src-bgd spectrum.
# We remove poor quality channels from the copy so we
# can make a plot which is consistent with the rgsflucer ouptut.
my $tmpSrcSpecFile = newFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, order => $order
, rgsorder => $order
, content => 'PLOTTABLE SOURCE SPECTRUM'
);
copyFile( source => $srcSpecFile , destination => $tmpSrcSpecFile );
my $data = readFITSTable(
file => $tmpSrcSpecFile
, extension => 'SPECTRUM'
, colname => [qw( COUNTS QUALITY )]
);
for my $i ( 0 .. $data->{-nrows} -1 )
{
$data->{COUNTS}[$i] = 0 if ($data->{QUALITY}[$i]);
if ( $order == 1 ){
# Total counts in plottable spectrum (quality filtered)
$total_spec_counts_1st_order += $data->{COUNTS}[$i];
# maximum effective integration time
$exposure = readFITSKeyword(
file => $tmpSrcSpecFile
, extension => 'SPECTRUM'
, keyword => 'EXPOSURE'
) or return exception();;
}
}
delete $data->{QUALITY};
delete $data->{-column}{QUALITY};
writeFITSColumn(
file => $tmpSrcSpecFile
, extension => 'SPECTRUM'
, column => $data
);
# Add to list of file for plotting
push( @spectra, $tmpSrcSpecFile );
my $totSpecFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, order => $order
, rgsorder => $order
, content => "RGS FITS SOURCE+BACKGROUND SPECTRUM"
);
doCommand(
'rgsspectrum', evlist => $evFile
, srclist => $srcList
, source => $srcnum
, rebin => 1
, order => $order
, edgechannels => 2
, withspectrum => 'Y'
, spectrumset => $totSpecFile
, bkgcorrect => 'N'
, withbkgset => 'N'
, spectrumbinning => 'lambda'
)
or return exception();
# Add to list of file for plotting
push( @totspectra, $totSpecFile );
doCommand( 'addattribute'
,set => [ $srcSpecFile , $totSpecFile ]
,attributename => ["SPECTRUM:SRCLABEL" , "SPECTRUM:SRC_RA" , "SPECTRUM:SRC_DEC"]
,attributetype => [qw(string real real)]
,attributecomment => [
'Source label'
,'Source RA'
,'Source DEC'
]
,stringvalue => $src_param{$srcnum}{LABEL}
,realvalue => [ $src_param{$srcnum}{RA} , $src_param{$srcnum}{DEC} ]
) or exception();
# Response matrix for order 1
if ( $CreateRgsResponseMatrix && ( ( $order == 1) || ( $order ==2 ) ) )
{
my $ResponseMatrix = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, order => $order
, rgsorder => $order
, content => "RGS RESPONSE MATRIX"
);
# Apply rgsrmfgen
doCommand(
'rgsrmfgen', rmfset => $ResponseMatrix
# , srclist => $RgsSrcList
, srclist => $srcList
, source => $srcnum
, spectrumset => $srcSpecFile
, evlist => $evFile
, order => $order
)
or return exception();
doCommand( 'addattribute'
,set => "$ResponseMatrix:MATRIX"
,attributename => ['SRCLABEL' , 'SRC_RA' , 'SRC_DEC']
,attributetype => [qw(string real real)]
,attributecomment => [
'Source label'
,'Source RA'
,'Source DEC'
]
,stringvalue => $src_param{$srcnum}{LABEL}
,realvalue => [ $src_param{$srcnum}{RA} , $src_param{$srcnum}{DEC} ]
) or exception();
push( @matrices, $ResponseMatrix );
}
}
# Plot the spectra
my $psFile = newFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, content => 'RGS SOURCE SPECTRUM'
, format => 'PS'
);
my $pdfFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, content => 'RGS SOURCE SPECTRUM'
, format => 'PDF'
);
doCommand(
'rgsspecplot', spectrumsets => [@spectra]
, sourcelistset => $srcList
, sourceid => $srcnum
, device => '/VPS'
, plotfile => $psFile
#, rebin => 'Y'
, rebin => 'N'
, mincounts => 20
)
or return exception();
# Convert to PostScript
PStoPDF(
source => $psFile
, destination => $pdfFile
)
or return exception();
## Begin rgslccorr testing
### # Set default time binning
### my $timebinsize = 1000;
# Set time binning:
info("RGS spectrum 1st order. Total counts : $total_spec_counts_1st_order");
info("RGS spectrum 1st order. Exposure time : $exposure");
my $timebinsize;
if ( $total_spec_counts_1st_order <= 0 ){
$timebinsize = 1000;
} else {
$timebinsize = 20.0 * ($exposure / $total_spec_counts_1st_order); # At least 20 counts per bin
$timebinsize = floor($timebinsize); # timebinsize must be integer
# 20 < $timebinsize < 1000
if ($timebinsize < 20){ $timebinsize = 20; }
if ($timebinsize > 1000){ $timebinsize = 1000; }
}
info("RGS spectrum 1st order. Time binning : $timebinsize");
# New file for RGS source light curve
my $bkglcname = '';
unless ($mode eq 'SpectoscopySmallWindow') {
$bkglcname = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, content => 'RGS SOURCE BACKGROUND TIMESERIES'
);
};
my $srclcname = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, content => 'RGS SOURCE TIMESERIES'
);
if ($mode eq 'SpectroscopySmallWindow') {
doCommand(
'rgslccorr', evlist => $evFile
, srclist => $srcList
, sourceid => $srcnum
, timebinsize => $timebinsize
, outputsrcfilename => $srclcname
, withbkgsubtraction => 'N'
) or return exception();
} else {
doCommand(
'rgslccorr', evlist => $evFile
, srclist => $srcList
, sourceid => $srcnum
, timebinsize => $timebinsize
, outputsrcfilename => $srclcname
, withbkgsubtraction => 'Y'
, outputbkgfilename => $bkglcname
) or return exception();
};
# Plot RGS TimeSeries if BACKV > 0. No plot otherwise
my ($extension, $colName) = ("RATE", "BACKV");
if ( hasFITSColumn( file => $srclcname, extension => $extension, column => $colName ) ){
# Check if all BACKV > 0
if ( ! ModuleUtil::hasColumnNegValues( $srclcname , $extension , $colName ) ){
# Start RGS TimeSeries plot (pdf)
info('Start elcplot for RGS ');
my $gtiExposureset = findFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'Exposure Attitude GTI'
);
info('Found gti for exposure ');
my $rgsTsPsFile = newFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, content => 'Source timeseries plot'
, format => 'PS'
);
info('Created ps file for RGS timeseries plot');
doCommand(
'elcplot'
, set => $srclcname
, usegtiset => 'yes'
, useflareset => 'no'
, gtiset => $gtiExposureset
, offset => 'yes'
, offsetstyle => 'auto'
, forcestart => 'yes'
, binsize => 1
, plotdevice => '/vcps'
, plotfile => $rgsTsPsFile
) or return exception();
my $rgstsPdfFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, src_num => $srcnum
, content => 'RGS SOURCE TIMESERIES PLOT'
, format => 'PDF'
);
PStoPDF(
source => $rgsTsPsFile
, destination => $rgstsPdfFile
) or return exception();
} else {
info("Skipping elcplot. Column $colName has negative values in $srclcname");
}
} else {
info("Skipping elcplot. No column $colName in $srclcname");
}
## end rgslccorr testing
} #foreach $srcnum
# Make a background model for each spectral order:
my %bkgModel = ();
foreach my $order ( 1, 2 )
{
# my $ok = doCommand('rgsbkgmodelTest'
# , table => "$evFile:EVENTS"
# );
# if ($ok && ($mode ne 'SpectroscopySmallWindow') )
my $rgsbkgTestCmd = "rgsbkgmodelTest table=$evFile:EVENTS -w 10";
my $rgsbkgTestCmdOut = Exec::run($rgsbkgTestCmd);
if ( (grep { /Background light curve cannot be computed/ } @{$rgsbkgTestCmdOut}) or ($mode eq 'SpectroscopySmallWindow') )
{
info("The data fail rgsbkgmodelTest or are SpectroscopySmallWindow. Skipping rgsbkgmodel.");
}
else
{
$bkgModel{$order} = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
, order => $order
, content => 'RGS BACKGROUND MODEL'
);
doCommand(
'rgsbkgmodel'
, table => "$evFile:EVENTS"
, order => $order
, backgroundmodelset => $bkgModel{$order}
, pdistincl => 95
)
or return exception();
}
}
# Extract the whole-field spectrum:
my $wholeFieldSrcFile = findFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
, content => 'RGS src list with whole field'
);
# my $inst = uc thisInstrument;
my $makeWholeFieldSpectrum = 0; # default
my $idxOfWhFS = 0; # default
if ($wholeFieldSrcFile)
{
my $srcdata = readFITSTable(
file => $wholeFieldSrcFile
, extension => 'SRCLIST'
, colname => [qw( INDEX LABEL )]
);
# Read source number and selection flag from RGS source list
my $idxCol = readFITSColumn( file => $wholeFieldSrcFile
, extension => 'SRCLIST'
, column => 'INDEX'
);
my $labelCol = readFITSColumn( file => $wholeFieldSrcFile
, extension => 'SRCLIST'
, column => 'LABEL'
);
for ( my $i=0 ; $i < $srcdata->{-nrows} ; $i++ )
{
next if( $srcdata->{LABEL}[$i] ne 'WHOLE_FIELD');
# Check that a region file exists for this source before selecting it.
my $idx = $srcdata->{INDEX}[$i];
unless ( hasFITSExtension( file => $wholeFieldSrcFile
, extension => "${inst}_SRC${idx}_SPATIAL"
)) {
info(
"Cannot find region files for RGS source $idx. No whole-field spectral products will be produced."
);
last;
}
$makeWholeFieldSpectrum = 1;
$idxOfWhFS = $idx;
}
}
if ($makeWholeFieldSpectrum)
{
my $spregion = "$wholeFieldSrcFile:${inst}_SRC${idxOfWhFS}_SPATIAL";
my @spectra;
my @totspectra;
# Create background+source spectrum for each order:
foreach my $order ( 1, 2 )
{
my $totSpecFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
# , src_num => 0
, order => $order
, rgsorder => $order
, content => "RGS WHOLE FIELD SRC+BKG SPECTRUM"
);
doCommand(
'rgsspectrum', evlist => $evFile
, srclist => $wholeFieldSrcFile
, source => $idxOfWhFS
, rebin => 1
, order => $order
, edgechannels => 2
, withspectrum => 'Y'
, spectrumset => $totSpecFile
, bkgcorrect => 'N'
, withbkgset => 'N'
, spectrumbinning => 'lambda'
)
or return exception();
# Add to list of file for plotting
push( @totspectra, $totSpecFile );
# Response matrix for order 1
if ( $CreateRgsResponseMatrix && ( ( $order == 1 ) || ( $order == 2 ) ) )
{
my $ResponseMatrix = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
# , src_num => 0
, order => $order
, rgsorder => $order
, content => "RGS WHOLE FIELD RESPONSE MATRIX"
);
# Apply rgsrmfgen
doCommand(
'rgsrmfgen', rmfset => $ResponseMatrix
, srclist => $wholeFieldSrcFile
, source => $idxOfWhFS
# , spectrumset => $srcSpecFile
, spectrumset => $totSpecFile
, evlist => $evFile
, order => $order
)
or return exception();
push( @matrices, $ResponseMatrix );
}
}
# Plot the total spectra
my $psFile = newFile(
class => 'intermediate'
, instrument => thisInstrument
, exp_id => $exp_id
# , src_num => 0
, content => 'whole field spectrum'
, format => 'PS'
);
my $pdfFile = newFile(
class => 'product'
, instrument => thisInstrument
, exp_id => $exp_id
# , src_num => 0
, content => 'RGS WHOLE FIELD SRC+BKG SPECTRUM'
, format => 'PDF'
);
doCommand(
'rgsspecplot', spectrumsets => [@totspectra]
, sourcelistset => $wholeFieldSrcFile
, sourceid => $idxOfWhFS
, device => '/VPS'
, plotfile => $psFile
#, rebin => 'Y'
, rebin => 'N'
, mincounts => 20
)
or return exception();
# Convert to PostScript
PStoPDF(
source => $psFile
, destination => $pdfFile
)
or return exception();
} #foreach $srcnum
return success();
}
1;