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.35'; $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 = '2018-04-19'; # ChangeLog # ========= # # version 2.35 - 2018-04-19 JVP # ------------ # # + RGS products processing is stopped if RGS presents the frames jump issue # ( number of entries mitmatch == STDGTI0[1-9] != EXPOSU0[1-9] ) # Info included in the Screening Report # # 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 ); #----------------------------------------------------------------------------------------------------- # Skip producing RGS Products if the evFile presents the known issue in the # RGS data == STDGTI0[1-9] != EXPOSU0[1-9] # my $rgs_issue_chk = 0; foreach my $i ( 1 .. 9){ my $gtiName = "STDGTI0" . $i; my $exposuName = "EXPOSU0" . $i ; if ( hasFITSExtension( file => $evFile0, extension => $gtiName ) && hasFITSExtension( file => $evFile0, extension => $exposuName )){ my $nrowsGTI = numberFITSRows( file => $evFile0, extension => $gtiName ); my $nrowsExposu = numberFITSRows( file => $evFile0, extension => $exposuName ); if ($nrowsGTI != $nrowsExposu){ my $instr = thisInstrument; info("RGS jump issue in $instr: $evFile0 constains difference number of rows in $gtiName and $exposuName extensions. This should not happen. Skipping further processing."); my $rgsIssueFile = newFile( class => 'intermediate' , instrument => thisInstrument , exp_id => $exp_id , content => 'RGS jumps issue' , format => 'ASCII' ); # Write something in the Screening Report # - Write checking files under intermediate/ dir my @lines = (); push (@lines, "RGS jump issue in $instr, Exp. $exp_id . Entries mitmatch in $gtiName and $exposuName\n"); writeASCIIFile( name => $rgsIssueFile , text => \@lines , append => 1 ); # This info (intermediate/*RGS_jumps_issue*) will be then included in the Screening Report -> script 'prep_screen_pm' calls 'mk_rgsspec_ps' $rgs_issue_chk = 1; # To stop RGS products processing last;} } } # Return to sucess() (no error) and do not create further products if the RGS issue is present return success() if ( $rgs_issue_chk ); #----------------------------------------------------------------------------------------------------- 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;