package ExpDetectPrep; 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 Regexp::Common qw(boolean); use Astro::FITS::CFITSIO qw(TRUE FALSE); # Declare identity, version, author, date, etc. $name = __PACKAGE__; $VERSION = '0.26'; $version = $VERSION; $author = 'Richard West,Dean Hinshaw,Duncan John Fyfe,Richard Saxton,Ian Stewart,Duncan Law-Green,Ed Chapin'; $date = '2013-09-05'; # # ChangeLog # ========= # # version 0.26 - 2013-09-05 # ------------------------- # # + Fix output TLMIN/MAX values after merge to preserve values calculated in EPICImageSize # # version 0.25 - 2012-12-12 # ------------------------- # # + Fix image event file merging for cases with 3+ exposures # # version 0.24 - 2010-09-07 # ------------------------- # # + Disable PN small window mode for source detection again (SRR request, see Mantis #17) # # version 0.23 - 2009-10-15 # ------------ # # + Use PN small window mode for source detection (SRR request, Mantis #17) # # version 0.22 - 2009-09-15 # ------------ # # + Do not use EPIC Open filter mode for source detection # # # version 0.21 - 2006-08-25 # ------------ # # + Do not use MOS mode PrimePartialW6 for source detection # # version 0.20 - 2006-08-16 # ------------ # # + DEBUGGING: Copy merged image after ONTIME update. Something is overwriting change. # # version 0.19 - 2006-08-16 # ------------ # # + Added test for PN low gain mode. Skip such exposures in source detection. # # version 0.18 - 2006-08-11 # ------------ # # + Old EXPOSURE update code was overwriting corrected values. # + Now select exposures based on EXPOSURE (good time) rather than exposure duration. # # version 0.17 - 2006-07-12 # ------------ # # + Merged image ONTIME correction changed to use addattribute. This will ensure # correction is applied by executable logs. # + XID merged exposure map was made as weighted addition of N bands. This # leaves the exposure map values too big. Use farith to divide weighted map by N to # get a better value. # # version 0.16 - 2006-05-18 # ------------ # # + When more than one exposure is selected for source detection # gtimerge the selected images and write the new ONTIME into the # merged image. This is so the proper ONTIME can be propogated to # the sources. # # version 0.15 - 2006-03-22 # ------------ # # + Corrected keepSafe parameters # # version 0.14 - 2006-02-21 # ------------ # # + Images supplied to the band 8 imweightadd merges are # written to a file to be used with emosaic in ImageMerge module. # Using the imweightadd output looses the indidual image keywords. # # version 0.13 - 2006-02-07 DJF # ------------ # # + Added test to make sure files exist before passing to imweightadd # + Fixed bad use of info() and foo(). This depends on info returning true. This is not guarenteed. # # version 0.12 - 2005-11-14 IMS # ------------ # # + Incorporated DJF's shortening of some intermediate file names. # # version 0.11 - 2005-11-04 IMS # ------------ # # + Slightly modified it to fix a bug which prevented band 9 images/exp maps from being made. # # version 0.10 - 2005-11-04 IMS # ------------ # # + Exposures for which flare screening was not performed are no longer considered for merging. # # version 0.09 - 2005-11-02 IMS # ------------ # # + Non-vig exposure maps (both those read and the newly created merged version) now have shorter content strings (thus shorter names). # + Merged exposure map now also made for band 8 (necessary for ImageMerge). # + Merged band 8 image is now made in the same way as the (new) merged band 8 exp map. # # version 0.08 - 2005-10-31 IMS # ------------ # # + XID images and exposure maps are now also being made (by adding merged images and exp maps for bands 2, 3 and 4). # # version 0.07 - 2005-10-28 IMS # ------------ # # + Altered things so that the band-8 imweightadd is not called if there are no input images. (I hope.) # # version 0.06 - 2005-10-21 RGW # ------------ # # + Make sure we look for FITS product images not intermediate files when creating the band 8 merged image # # version 0.05 - 2005-10-14 IMS # ------------ # # + Merged image (as intermediate file) now made for band 8 as requested by ACS. # + Added copyFile to the list of used ModuleResources subs. # # version 0.04 - 2005-09-21 IMS # ------------ # # + Slight bug fixing, so it will run. # # version 0.03 - 2005-09-20 IMS # ------------ # # + Restored imweightadd parameters to 0.01 format. There's method in my madness. # + Constituent exposure ids are now written as kwds to merged datasets. # # version 0.02 - 2005-08-15 RGW # ------------ # # + hacked argument to imweightadd (withweights -> calculateweights) # + changed circumstances under which the module is ignored (only ignored if # all image creation modules are ignored) # # version 0.01 - 2005-06-12 IMS # ------------ # # + First draft. # Declare list of instruments this module is interested in @instList = qw(emos1 emos2 epn); #my $energyBands = [qw(1 2 3 4 5)]; my $xidBands = {'epn' => [qw(2 3 4)], 'emos1' => [qw(2 3 4)], 'emos2' => [qw(2 3 4)] }; ### This list of xid bands ought to be obtained from a subroutine call because it is also needed in ExpDetect. IMS 2005-10-31. # Number of streams sub numberOfStreams { 1; } # sub evaluateRules { my $inst = thisInstrument; if ( $inst eq 'epn' ) { # PN Ignore condition return ignore() if ( allIgnored( module => 'MakePNImage' , instrument => thisInstrument ) ); # PN Start conditions start() if ( allComplete( module => 'MakePNImage' , instrument => thisInstrument ) ); } else { # MOS Ignore condition return ignore() if ( allIgnored( module => 'MakeMOSImage' , instrument => thisInstrument ) ); # MOS Start conditions start() if ( allComplete( module => 'MakeMOSImage' , instrument => thisInstrument ) ); } } # sub performAction { info("Module version number: $version"); my @exp_id = listExposures( instrument => thisInstrument ); # Set which bands to use for source searching my @bandList = ( 1, 2, 3, 4, 5, 8 ); # If more than 1 exposure, possibly with different filters, # find which filter has the longest total exposure: my %exposureLengths; my (%filters, %modes); foreach my $exp_id (@exp_id) { my @opt; my $mode = getExposureProperty( instrument => thisInstrument , exp_id => $exp_id , name => 'mode' ); info("DEBUG: Exposure $exp_id Mode $mode"); ## Skip EPIC PN|MOS Open filter for source detection if ( thisInstrument =~ /epn|mos/ ) { my $filter = getExposureProperty( instrument=> thisInstrument , exp_id => $exp_id , name => 'filter' ); if ( $filter =~ /Open/ ) { info("EPIC Open filter should not be used for source detection"); next; } } if ( thisInstrument eq 'epn' ) { # EC (14 NOV 2013): re-enable at PR's request -- comment out the code that skips # Skip small window mode. # Re-enabled at SRR's request (SEE: Mantis #17) #if ( $mode =~ /Small/ ) #{ # info("Exposure $exp_id: Will not perform source detection for PN mode $mode"); # next; #} # Skip low gain mode my $evFile=findFile( class => 'product' , instrument => thisInstrument , exp_id => $exp_id , content => 'EPIC PN imaging mode event list' ); if ( $evFile && ModuleUtil::pnLowGainWarn($evFile) ) { info("Exposure $exp_id: EPN exposure in low gain mode - won't perform source detection on it"); next; } } elsif ( thisInstrument =~ /emos/ && $mode =~ /PrimePartialW6/i ) { info("Exposure $exp_id: Will not perform source detection for MOS mode $mode"); next; } my $rawim = findFile( class => 'product' , instrument => thisInstrument , exp_id => $exp_id , band => 1 , content => 'EPIC image' , 'format' => 'FITS' ); if ( !$rawim ) { info("Exposure $exp_id: Can't find a band 1 image for this exposure - won't perform source detection on it"); next; } my $flareScreened; @opt = ( file => $rawim , extension => 'PRIMARY' , keyword => 'BKGDSCRN'); if ( hasFITSKeyword( @opt ) ) { $flareScreened = readFITSKeyword( @opt ); $flareScreened = ($flareScreened =~ /$RE{boolean}{true}/i) ? 1 : 0 ; } if ( !$flareScreened ) { info("Images were not flare screened - source detection not performed"); next; } my $filter = getExposureProperty( instrument => thisInstrument , exp_id => $exp_id , name => 'filter' ); unless ( $filter =~ /(Thin|Medium|Thick|Open)/ ) { info("Could not match filter for energy conversion factor ($filter)"); next; # return exception(); } $filter = $1; $filters{$exp_id} = $filter; $modes{$exp_id} = $mode; my $exposureLength = 0; @opt = ( file => $rawim , extension => 'PRIMARY' , keyword => 'EXPOSURE' ) ; if ( hasFITSKeyword( @opt ) ) { $exposureLengths{$filter}{$mode} += readFITSKeyword( @opt ); } } # end loop over exposure ids return success() if (keys(%exposureLengths)<=0); my $maxExposureLength = 0; my ($chosenFilter, $chosenMode); foreach my $filter (keys(%exposureLengths)) { foreach my $mode (keys(%{$exposureLengths{$filter}})) { if ($exposureLengths{$filter}{$mode} > $maxExposureLength) { $maxExposureLength = $exposureLengths{$filter}{$mode}; $chosenFilter = $filter; $chosenMode = $mode; } } } info("Chosen filter $chosenFilter"); info("Chosen mode $chosenMode"); my $tempset = newFile( class => 'intermediate' , content => 'Generic temporary dataset' ); foreach my $band (@bandList) { my (@rawImList, @expImList,@expidList); foreach my $exp_id (@exp_id) { next unless (($filters{$exp_id} eq $chosenFilter) && ( $modes{$exp_id} eq $chosenMode ) ); my $rawim = findFile( class => 'product' , instrument => thisInstrument , exp_id => $exp_id , band => $band , content => 'EPIC image' , 'format' => 'FITS' ); my $expim = findFile( class => 'product' , instrument => thisInstrument , exp_id => $exp_id , band => $band , content => 'EPIC exposure map' , 'format' => 'FITS' ); next unless ( $rawim && fileExists( file => $rawim ) && $expim && fileExists( file => $expim ) ); push( @rawImList, $rawim ); push( @expImList, $expim ); push( @expidList, $exp_id ); } # end loop over $exp_id unless ( @rawImList ) { info("No images found"); next; } foreach my $exp_id (@expidList) { setExposureProperty( instrument => thisInstrument , exp_id => $exp_id , name => 'srcdet' , value => 0+TRUE ); } info("Per camera merging of images for ",thisInstrument); my $mergedRawim = newFile( class => 'intermediate' , instrument => thisInstrument , band => $band , content => 'merged image' ); doCommand( 'imweightadd' , imagesets => [@rawImList] , withweights => 'no' , outimageset => $mergedRawim , tempset => $tempset ) or return exception(); # When merging more than one exposure calculate new ONTIME and EXPOSURE values # and write this into the merged image. if (@rawImList > 1) { info("Calculating corrected ONTIME for merged images"); } my $totalExposure; my %ontime = ModuleUtil::ontimeFix( $mergedRawim , @rawImList ); if (my @ontime = sort( grep /ONTIME/,(keys %ontime)) ) { my $num = scalar(@ontime); my $numnum = $num+1; my @attcomment = ('Sum of GTIs for a particular CCD') x $num; push @attcomment , 'max of ONTIMEnn values' ; my @wattcomment = ('Y') x $numnum; my @atttype = ('real') x $numnum; my @attname = (@ontime , 'EXPOSURE'); my @realval = (@ontime{@ontime} , $ontime{EXPOSURE}); ### DJF 2006-08-18 ### For some reason unknown to me the first addattribute below is recorded in the $mergedRawim ### HISTORY but the ONTIME and EXPOSURE updates are lost.the ONTIME and EXPOSURE updates are done correctly. doCommand('addattribute' , set => $mergedRawim ,attributename => \@attname ,attributetype => \@atttype ,attributecomment => \@attcomment ,realvalue => \@realval ) or return exception(); doCommand('addattribute' , set => $mergedRawim ,attributename => \@attname ,attributetype => \@atttype ,attributecomment => \@attcomment ,realvalue => \@realval ) or return exception(); $totalExposure = $ontime{EXPOSURE}; } # Write the filenames merged to a file so we can use them # again with emosaic. my @rawlistdtl = ( class => 'intermediate' , content => 'raw image list' , instrument => thisInstrument , band => $band , format => 'ASCII' ); unless ( ModuleUtil::keepSafe ( file_identifier => \@rawlistdtl , text => \@rawImList ) ) { info("Unable to store raw image list for later use."); return exception(); } foreach my $i (0 .. $#expidList) { my $exp_id = $expidList[$i]; my $kwd_str = 'EXPID'.sprintf("%2.2d", $i+1); doCommand( 'addattribute', set => $mergedRawim.":PRIMARY" , attributename => $kwd_str , attributetype => 'string' , stringvalue => $exp_id ) or return exception(); } my $mergedExpim = newFile( class => 'intermediate' , instrument => thisInstrument , band => $band , content => 'merged exposure map' ); doCommand( 'imweightadd' , imagesets => [@expImList] , withweights => 'no' , outimageset => $mergedExpim , tempset => $tempset ) or return exception(); # Write the filenames merged to a file so we can use them # again with emosaic. my @explistdtl = ( class => 'intermediate' , content => 'exposure map list' , instrument => thisInstrument , band => $band , format => 'ASCII' ); unless ( ModuleUtil::keepSafe ( file_identifier => \@explistdtl , text => \@expImList ) ) { info("Unable to store exposure map list for later use."); return exception(); } doCommand( 'addattribute', set => $mergedExpim.":PRIMARY" , attributename => 'EXPOSURE' , attributetype => 'real' , realvalue => $totalExposure ) or return exception(); foreach my $i (0 .. $#expidList) { my $exp_id = $expidList[$i]; my $kwd_str = sprintf('EXPID%02d' , $i+1); doCommand( 'addattribute', set => $mergedExpim.":PRIMARY" , attributename => $kwd_str , attributetype => 'string' , stringvalue => $exp_id ) or return exception(); } # next if ($band eq '8'); # band 8 is just for rawimage and expmap. my (@flatExpImList, @imEventList); foreach my $exp_id (@exp_id) { next unless ($filters{$exp_id} eq $chosenFilter && $modes{$exp_id} eq $chosenMode ); my $flatExpIm = findFile( class => 'intermediate' , instrument => thisInstrument , exp_id => $exp_id , band => $band , content => 'non-vig exposure map' ); my $filteredList = findFile( class => 'intermediate' , instrument => thisInstrument , exp_id => $exp_id , band => $band , content => 'image event list' ); next unless $flatExpIm && $filteredList; push( @flatExpImList, $flatExpIm ); push( @imEventList, $filteredList ); } # end loop over $exp_id unless (@flatExpImList) { info("No flat expmaps found"); next; } my $mergedFlatExpim = newFile( class => 'intermediate' , instrument => thisInstrument , band => $band , content => 'merged non-vig exp map' ); doCommand( 'imweightadd' , imagesets => [@flatExpImList] , withweights => 'no' , outimageset => $mergedFlatExpim , tempset => $tempset ) or return exception(); # Write the filenames merged to a file so we can use them # again with emosaic. doCommand( 'addattribute', set => $mergedFlatExpim.":PRIMARY" , attributename => 'EXPOSURE' , attributetype => 'real' , realvalue => $totalExposure ) or return exception(); foreach my $i (0 .. $#expidList) { my $exp_id = $expidList[$i]; my $kwd_str = 'EXPID'.sprintf("%2.2d", $i+1); doCommand( 'addattribute', set => $mergedFlatExpim.":PRIMARY" , attributename => $kwd_str , attributetype => 'string' , stringvalue => $exp_id ) or return exception(); } my $mergedFilteredList = newFile( class => 'intermediate' , instrument => thisInstrument , band => $band , content => 'image merged event list' ); my $imEventListSize = scalar(@imEventList); info("DEBUG: + + + Number of image event lists: $imEventListSize"); # Here we will merge together multiple events lists if needed # (Branch>1). merge messes up the TLMIN/MAX keywords that are # set in EPICImageSize to handle mosaics. If the input events # lists have identical TL[MIN/MAX][6/7] values they will be # re-written upon completion to the merged output list as a # hack. my $fixtlm=1; my ($tlmin6,$tlmax6,$tlmin7,$tlmax7); if (scalar(@imEventList)>0) { $tlmin6 = int(readFITSKeyword(file => $imEventList[0], extension => 'EVENTS', keyword => 'TLMIN6')); $tlmax6 = int(readFITSKeyword(file => $imEventList[0], extension => 'EVENTS', keyword => 'TLMAX6')); $tlmin7 = int(readFITSKeyword(file => $imEventList[0], extension => 'EVENTS', keyword => 'TLMIN7')); $tlmax7 = int(readFITSKeyword(file => $imEventList[0], extension => 'EVENTS', keyword => 'TLMAX7')); } if (scalar(@imEventList)==1) { info("DEBUG: Branch 1"); $fixtlm=0; copyFile( source => $imEventList[0] , destination => $mergedFilteredList ); } elsif(scalar(@imEventList)==2) { info("DEBUG: Branch 2"); if( ($tlmin6 != int(readFITSKeyword(file => $imEventList[1], extension => 'EVENTS', keyword => 'TLMIN6'))) || ($tlmax6 != int(readFITSKeyword(file => $imEventList[1], extension => 'EVENTS', keyword => 'TLMAX6'))) || ($tlmin7 != int(readFITSKeyword(file => $imEventList[1], extension => 'EVENTS', keyword => 'TLMIN7'))) || ($tlmax7 != int(readFITSKeyword(file => $imEventList[1], extension => 'EVENTS', keyword => 'TLMAX7'))) ) { $fixtlm=0; info("Interesting: $imEventList[1] has different TL[MIN/MAX][6/7] keywords than $imEventList[0]. Won't fix merge output."); } doCommand( 'merge' , set1 => $imEventList[0] , set2 => $imEventList[1] , outset => $mergedFilteredList ) or return exception(); } else { info("DEBUG: Branch 3+"); my $prevFile=shift(@imEventList); my $idx=0; # foreach my $file (splice(@imEventList, 2)) foreach my $file (@imEventList) { info("DEBUG: > > > Step $idx"); my $nextFile=newFile(class => 'intermediate', content => 'Intermediate merged event list', src_num => $idx ); $idx++; if( ($tlmin6 != int(readFITSKeyword(file => $file, extension => 'EVENTS', keyword => 'TLMIN6'))) || ($tlmax6 != int(readFITSKeyword(file => $file, extension => 'EVENTS', keyword => 'TLMAX6'))) || ($tlmin7 != int(readFITSKeyword(file => $file, extension => 'EVENTS', keyword => 'TLMIN7'))) || ($tlmax7 != int(readFITSKeyword(file => $file, extension => 'EVENTS', keyword => 'TLMAX7'))) ) { $fixtlm=0; info("Interesting: $file has different TL[MIN/MAX][6/7] keywords than $imEventList[0]. Won't fix merge output."); } doCommand('merge', ,set1 => $prevFile ,set2 => $file ,outset => $nextFile ) or return exception(); $prevFile=$nextFile; } copyFile(source => $prevFile, destination => $mergedFilteredList); } if ($fixtlm) { info("Merge may have broken TL[MIN/MAX}[6/7] keywords. Hacking $mergedFilteredList: TLMIN6=$tlmin6 TLMAX6=$tlmax6 TLMIN7=$tlmin7 TLMAX7=$tlmax7"); FITSUtils::writeKeyword($mergedFilteredList, "EVENTS", "TLMIN6" => {value => $tlmin6}); FITSUtils::writeKeyword($mergedFilteredList, "EVENTS", "TLMAX6" => {value => $tlmax6}); FITSUtils::writeKeyword($mergedFilteredList, "EVENTS", "TLMIN7" => {value => $tlmin7}); FITSUtils::writeKeyword($mergedFilteredList, "EVENTS", "TLMAX7" => {value => $tlmax7}); } foreach my $i (0 .. $#expidList) { my $exp_id = $expidList[$i]; my $kwd_str = 'EXPID'.sprintf("%2.2d", $i+1); doCommand( 'addattribute', set => $mergedFilteredList , attributename => $kwd_str , attributetype => 'string' , stringvalue => $exp_id ) or return exception(); } } # end loop over $band # Merged image and exposure map for band 9 (xid band). These images should be used in the second call to emldetect in ExpDetect. my @rawImList=(); my @expImList=(); my $inst = thisInstrument; my @xidBandsThisInst = @{$xidBands->{$inst}}; foreach my $band (@xidBandsThisInst) { my $rawim = findFile( class => 'intermediate' , instrument => thisInstrument , band => $band , content => 'merged image' ); my $expim = findFile( class => 'intermediate' , instrument => thisInstrument , band => $band , content => 'merged exposure map' ); unless ( $rawim && fileExists( file => $rawim ) && $expim && fileExists( file => $expim ) ) { info("Band $band image not found"); next; } push( @rawImList, $rawim ); push( @expImList, $expim ); } # end loop over $band unless (@rawImList == @xidBandsThisInst) { info("\nNo or missing XID-band constituent images."); return success(); } $tempset = newFile( class => 'intermediate' , content => 'Generic temporary dataset' ); my $mergedRawim = newFile( class => 'intermediate' , instrument => thisInstrument , band => 9 , content => 'merged image' ); doCommand( 'imweightadd' , imagesets => [@rawImList] , withweights => 'no' , outimageset => $mergedRawim , tempset => $tempset ) or return exception(); my $addedExpim = newFile( class => 'intermediate' , instrument => thisInstrument , band => 9 , content => 'added exposure map' ); my $mergedExpim = newFile( class => 'intermediate' , instrument => thisInstrument , band => 9 , content => 'merged exposure map' ); doCommand( 'imweightadd' , imagesets => [@expImList] , withweights => 'no' , outimageset => $addedExpim , tempset => $tempset ) or return exception(); my $div = @expImList; doCommand('farith' ,infil1 => $addedExpim ,infil2 => $div ,outfil => $mergedExpim ,ops => 'DIV' ) or return exception(); # Raise success flag return success(); } 1;