package Subs::UvotImages; ############################################################################## # # DESCRIPTION: Extract transform and plot images # # HISTORY: # HISTORY: $Log: UvotImages.pm,v $ # HISTORY: Revision 1.89 2011/01/20 18:39:26 apsop # HISTORY: Added new filter usage when calling sextractor # HISTORY: # HISTORY: Revision 1.86 2009/12/18 14:37:35 apsop # HISTORY: This release contains several new items: # HISTORY: 1. call to Subs::UvotGraspCorr for aspect correction # HISTORY: 2. call to Subs::UvotProduct, which generates gifs light curves # HISTORY: 3. call to Subs::XrtGrbLc, which generates gifs light curves for GRB based on XRT data # HISTORY: # HISTORY: Revision 1.85 2008/07/28 15:38:28 apsop # HISTORY: Add max.rate=1000 to starid in uvotskycorr. Always use SHIFTADD method in uvotexpmap. # HISTORY: # HISTORY: Revision 1.84 2008/06/23 16:38:17 apsop # HISTORY: Check for existence of files before applying aspect correction. # HISTORY: # HISTORY: Revision 1.83 2008/05/16 15:01:17 apsop # HISTORY: Round real image values. Use SHIFTADD method in uvotexpmap where we can. Fix bug in setting of ATTSTATU keyword. Run uvotexpmap on non-aspect-corrected images. # HISTORY: # HISTORY: Revision 1.82 2007/12/20 19:43:20 apsop # HISTORY: Fix bug in setting of exclude param in uvotimsum. # HISTORY: # HISTORY: Revision 1.81 2007/12/18 19:59:02 apsop # HISTORY: Make rate images from count images and use to make plots. # HISTORY: # HISTORY: Revision 1.80 2007/11/08 17:10:13 apsop # HISTORY: Do not try and use unmask file with grism data. # HISTORY: # HISTORY: Revision 1.79 2007/09/28 17:14:48 apsop # HISTORY: Change writing of ATTFLAG keyword so that it is a string rather than an integer. # HISTORY: # HISTORY: Revision 1.78 2007/09/17 15:34:03 apsop # HISTORY: Fix bug in determining mask file filename. # HISTORY: # HISTORY: Revision 1.77 2007/09/12 18:03:30 apsop # HISTORY: Protect against non-existent intermediate image some files. Also protect against missing catalogue housekeeping file. # HISTORY: # HISTORY: Revision 1.76 2007/09/11 18:08:38 apsop # HISTORY: Pass mask file from uvotexpmap to uvotimsum. Use att jump corrected attitude # HISTORY: file in uvot image processing. Find AT settling image and set OBS_MODE to # HISTORY: SETTLING. # HISTORY: # HISTORY: Revision 1.75 2007/06/28 20:37:56 apsop # HISTORY: Fix bugs in support for indexed image files. # HISTORY: # HISTORY: Revision 1.74 2007/04/18 22:42:44 apsop # HISTORY: Rework image processing to allow for indexed image files. # HISTORY: # HISTORY: Revision 1.73 2007/03/23 13:05:43 apsop # HISTORY: Turn off flat fielding in uvotimage. # HISTORY: # HISTORY: Revision 1.72 2007/03/15 21:29:20 apsop # HISTORY: Do mod8 and flat field correction for final processing. # HISTORY: # HISTORY: Revision 1.71 2007/02/08 21:34:51 apsop # HISTORY: Write ATTFLAG keyword to all extensions of corrected file. # HISTORY: # HISTORY: Revision 1.70 2007/02/02 17:19:11 apsop # HISTORY: Change parameters for uvotdetect. # HISTORY: # HISTORY: Revision 1.69 2007/02/01 15:06:44 apsop # HISTORY: Run uvotattcorr to produce corrected attitude file. Changes to the look of the uvot image plots. # HISTORY: # HISTORY: Revision 1.68 2006/08/21 23:06:12 apsop # HISTORY: Removed the code for figuring out which HDUs to exclude from the summed # HISTORY: images. As of Build 18 the HDUs without aspect corrections are excluded # HISTORY: by default (and a bug was preventing passing the list from working). # HISTORY: # HISTORY: Revision 1.67 2006/07/28 13:17:18 apsop # HISTORY: Don;t do any processing on BLOCKED images or if event files. # HISTORY: # HISTORY: Revision 1.66 2006/07/03 01:47:20 apsop # HISTORY: Fix bug in construction of gif plot titles. # HISTORY: # HISTORY: Revision 1.65 2006/06/28 19:09:41 apsop # HISTORY: Aesthetic changes for image gif files: no smoothing, black body color, # HISTORY: better titles. # HISTORY: # HISTORY: Revision 1.64 2006/05/10 15:18:29 apsop # HISTORY: Change filename for aspect correction output to new housekeeping file. # HISTORY: # HISTORY: Revision 1.63 2006/05/06 21:12:49 apsop # HISTORY: Put show command in ximage to track version. *_PNT keywords gotten from # HISTORY: ra,dec,roll parameters instead of burst_* parameters. # HISTORY: # HISTORY: Revision 1.62 2006/04/26 20:48:12 apsop # HISTORY: Replace partition parameter with catspec parameter. # HISTORY: # HISTORY: Revision 1.61 2006/03/13 17:45:59 apsop # HISTORY: Make error message about no data more explicit. # HISTORY: # HISTORY: Revision 1.60 2005/12/22 18:05:09 apsop # HISTORY: Add * for index number in Filename::get method, so that all unfiltered event and # HISTORY: raw image files are picked up. # HISTORY: # HISTORY: Revision 1.59 2005/12/16 19:27:44 apsop # HISTORY: Pass attitude to swiftxform when creating grism DET images to get # HISTORY: alternate SKY WCS keys. # HISTORY: # HISTORY: Revision 1.58 2005/12/15 22:37:54 apsop # HISTORY: Pass the expmap flag to uvotimsum to indicate whether or not exposure # HISTORY: maps are being summed. # HISTORY: # HISTORY: Revision 1.57 2005/11/15 13:52:17 apsop # HISTORY: Fix inexplicablly uncaught bug in in mispelling of unlink. # HISTORY: # HISTORY: Revision 1.56 2005/11/09 20:13:56 apsop # HISTORY: Get star catalogue location from sw0.par parameter starcatalog. # HISTORY: # HISTORY: Revision 1.55 2005/11/09 16:30:45 apsop # HISTORY: don;t try and do any image processing if there are no image or event files. # HISTORY: # HISTORY: Revision 1.54 2005/11/08 20:00:18 apsop # HISTORY: <Previous comment bogus>Clean up temporary input list files. # HISTORY: # HISTORY: Revision 1.53 2005/11/08 19:22:28 apsop # HISTORY: Populate the TIMELIST and DATALIST hashes. Used to be an SWCheckInput. # HISTORY: # HISTORY: Revision 1.52 2005/11/03 21:33:01 apsop # HISTORY: Modified a few parameters to work reasonably with Swift Build 16. # HISTORY: # HISTORY: Revision 1.51 2005/11/02 16:34:20 apsop # HISTORY: Use uvotimage to create Level 1 and 2 images. Use uvotskycorr to find # HISTORY: aspect corrections. Only include aspect corrected images in sums. # HISTORY: # HISTORY: Revision 1.50 2005/08/30 14:13:47 apsop # HISTORY: Copy date keywords from sky image prime header to exposure map. # HISTORY: # HISTORY: Revision 1.49 2005/06/01 13:55:27 apsop # HISTORY: Test for presence of event file before trying to use it. # HISTORY: # HISTORY: Revision 1.48 2005/04/29 15:47:24 apsop # HISTORY: Disable calling of uvotmodmap. # HISTORY: # HISTORY: Revision 1.47 2005/04/06 15:43:42 apsop # HISTORY: Change to using CALDB for cal parameters. # HISTORY: # HISTORY: Revision 1.46 2005/03/16 13:31:37 apsop # HISTORY: Add in FILTER keyword to image primary header. # HISTORY: # HISTORY: Revision 1.45 2005/03/15 18:58:39 apsop # HISTORY: Put basic and DATE* keywords into image file primary headers. Fix bug (?) # HISTORY: which was not using modmap output for processing. # HISTORY: # HISTORY: Revision 1.44 2005/02/08 18:21:32 apsop # HISTORY: Add new types for uvot filter specific images and exp maps. Replace # HISTORY: Filename::corresponding method with file name hack for now. Fix bug which was # HISTORY: preventing the production of raw images. # HISTORY: # HISTORY: Revision 1.43 2005/01/21 04:26:43 apsop # HISTORY: Change to use file names stored in cat file rather than contructing them. # HISTORY: # HISTORY: Revision 1.42 2005/01/12 17:29:23 apsop # HISTORY: Change weightfile param in uvotdetect to expfile param. # HISTORY: # HISTORY: Revision 1.41 2005/01/08 04:11:00 apsop # HISTORY: Change test for existence of event list to be more forgiving. # HISTORY: # HISTORY: Revision 1.40 2005/01/07 20:34:21 apsop # HISTORY: No longer necessary to use the CENTER method of swiftxform or a small # HISTORY: value of ncell for uvotmodmap. # HISTORY: # HISTORY: Revision 1.39 2004/12/03 13:38:11 apsop # HISTORY: Turn back on calling uvotmodmap. # HISTORY: # HISTORY: Revision 1.38 2004/10/13 01:38:31 apsop # HISTORY: Disable uvotmodmap for the time being, and fix bug in getting plot file names. # HISTORY: # HISTORY: Revision 1.37 2004/10/12 16:26:36 apsop # HISTORY: Filename changes and a few other tweaks to make module run in complete pipeline. # HISTORY: # HISTORY: Revision 1.36 2004/09/15 22:26:24 apsop # HISTORY: Remove the binning factor from uvot event list file names # HISTORY: # HISTORY: Revision 1.35 2004/09/01 14:42:21 apsop # HISTORY: Disable level 2/3 processing for now. # HISTORY: # HISTORY: Revision 1.34 2004/08/31 20:51:02 apsop # HISTORY: Implemented UVOT level 2 and 3 products. # HISTORY: # HISTORY: Revision 1.33 2004/08/18 17:27:58 apsop # HISTORY: Yet another fix for the timing keywords. # HISTORY: # HISTORY: Revision 1.32 2004/07/19 16:03:24 apsop # HISTORY: Add version number entry to comment field. # HISTORY: # HISTORY: Revision 1.31 2004/07/12 13:46:09 apsop # HISTORY: Fix problems with timing keywords # HISTORY: # HISTORY: Revision 1.30 2004/06/15 16:01:04 apsop # HISTORY: Better trapping of problem uvot exposures. # HISTORY: # HISTORY: Revision 1.29 2004/06/02 16:11:19 apsop # HISTORY: Change name of catalogue extension. # HISTORY: # HISTORY: Revision 1.28 2004/05/04 16:33:01 dah # HISTORY: Bug fixes for case where no initial image file exists. # HISTORY: # HISTORY: Revision 1.27 2004/04/30 16:16:40 dah # HISTORY: Use new FITSfile::keywords() method to speed up image extraction # HISTORY: # HISTORY: Revision 1.26 2004/04/16 20:21:18 dah # HISTORY: Begin using embedded history records # HISTORY: # HISTORY: Revision 1.25 2004/04/16 19:53:18 dah # HISTORY: Begin using embedded history records # HISTORY: # # VERSION: 0.0 # ############################################################################## use Subs::Sub; use Subs::Images; use Util::Xanadu; use Subs::UvotNames; @ISA = ("Subs::Images"); use strict; sub new { my $proto=shift; my $self=$proto->SUPER::new(); $self->{DESCRIPTION}="Extracting, merging and plotting images for UVOT"; return $self; } ################## # METHODS: ################## sub body { my $self=shift; my $log =$self->log(); my $filename=$self->filename(); my $catfile = $filename->get('hk', 'uvot', 'ct', '*'); unless( -f $catfile ){ $log->entry("No uvot catalogue file indicates no uvot data"); return; } $self->create_images(); $self->image_raw_to_det(); $self->create_exposure_maps(); $self->aspect_correct_sky_images(); $self->tag_round_images(); # create summed images/exposure maps $self->sum_images(); # detect sources and determine magnitudes $self->detect_sources(); $self->plot_images(); } # end of body method ############################################################################# # Extract images from filtered event data and raw images files ############################################################################# sub create_images { my $self = shift; my $log =$self->log(); my $filename=$self->filename(); my $procpar =$self->procpar(); my $jobpar =$self->jobpar(); $log->entry('Creating images from UVOT event and image mode data'); my $attfile = $filename->get('attitude', 's'); my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if(-f $uattfi){ $jumpatt = $uattfi; } $attfile = $jumpatt if -f $jumpatt; if (not $attfile) { $log->error(2, 'no attitude file available'); return; } my @eventFiles = $filename->get('unfiltered', 'uvot', '*', '*'); @eventFiles = grep !/[0-9]ubl/, @eventFiles; my @imageFiles = $filename->get('rawimage', 'uvot', '*', '*'); @imageFiles = grep !/[0-9]ubl/, @imageFiles; @imageFiles = grep !/[0-9]udi/, @imageFiles; unless( @eventFiles || @imageFiles ){ $log->error(1, 'No UVOT image or event files to process.'); return; } my $catfile = $filename->get('hk', 'uvot', 'ct', '*'); if (not -f $catfile) { $log->entry('No UVOT exposure catalogue to update'); }else { $log->entry('UVOT exposure catalogue update not implemented'); } my $prefix = 'Qz3'; # a unique string my $tool = Util::HEAdas->new('uvotimage')->is_script(1); $tool->params({ prefix => $prefix, attfile => $attfile, teldeffile => 'CALDB', alignfile => 'CALDB', ra => $jobpar->read('ra'), dec => $jobpar->read('dec'), roll => $jobpar->read('roll'), flatfield => 'no', mod8corr => 'no', # badpix => 'no', # catfile => $catfile, chatter => 5 }); $tool->params({mod8corr => 'yes'}) if $jobpar->{TIMELIST}->{final}; my $swobsid = $filename->sequence_specific; my @dataFiles = (@eventFiles, @imageFiles); while( @dataFiles ){ my $type = $dataFiles[0] =~ /\.evt/ ? 'unfiltered' : 'rawimage'; my ($mode, $index) = ($filename->parse($dataFiles[0], $type))[1,2]; $mode = substr($mode, 0, 2); my @infiles = ( $filename->get('rawimage', 'uvot', "$mode*", $index), $filename->get('unfiltered', 'uvot', "$mode*", $index) ); my $expr = '('. join('|',@infiles) .')'; @dataFiles = grep !/${expr}/, @dataFiles; $tool->params({infile => join(',', @infiles)})->run(); # move each file ${prefix}xxx to sw<obsid>xxx foreach my $name (glob($prefix . '*')) { my $xxx = substr($name, length($prefix)); my $fixed = $swobsid . $xxx; if($index!=0){ $index = sprintf('%02d', $index) if length($index)<2; $fixed =~ s/\./_${index}./; } if (-f $fixed) { unlink($fixed); } rename($name, $fixed); } } } # end create images method ############################################################################### # convert raw coordinate images to detector coordinates ############################################################################### sub image_raw_to_det { my $self=shift; my $log =$self->log(); my $filename=$self->filename(); my $procpar =$self->procpar(); my $jobpar =$self->jobpar(); $log->entry("Converting raw grism images to detector coordinates"); my $attfile = $filename->get('attitude', 's'); my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if (-f $uattfi) { $jumpatt = $uattfi; } $attfile = $jumpatt if -f $jumpatt; if (not $attfile) { $attfile = 'NONE'; } $log->entry("Processing raw uvot images"); my $ubadpix = Util::HEAdas->new('uvotbadpix'); my $umodmap = Util::HEAdas->new('uvotmodmap'); my $uflatfield = Util::HEAdas->new('uvotflatfield'); $uflatfield->params({flatfile => 'CALDB'}); $log->entry("setting {RA,DEC,PA}_PNT in all raw images"); my $ra = $jobpar->read("ra"); my $dec = $jobpar->read("dec"); my $roll = $jobpar->read("roll"); ######################################### # create the swiftxform tool ######################################### my $swiftxform = Util::HEAdas->new('swiftxform'); $swiftxform->params({ teldeffile => 'CALDB', to => 'DET', attfile => $attfile, ra => $jobpar->read('ra'), dec => $jobpar->read('dec'), roll => $jobpar->read('roll'), aberration => 'no', seed => $procpar->read('seed'), chatter => 4, clobber => 'yes', }) ->is_script(1); ###################### # loop over _GRISM_ files # since only they are converted to DET coordinates ###################### foreach my $rawFile ($filename->get('rawimage', 'uvot', 'g*', '*') ) { { my $imageFits = Util::FITSfile->new($rawFile); for my $i (1 .. $imageFits->nhdus - 1) { $imageFits->ext($i); $imageFits->keyword('RA_PNT', $ra); $imageFits->keyword('DEC_PNT', $dec); $imageFits->keyword('PA_PNT', $roll); } } my $badFile = $filename->corresponding('rawimage', 'badimage', $rawFile); my $corrFile = $filename->corresponding('rawimage', 'corrimage', $rawFile); my $mod8File = 'mod.tmp'; $log->entry("running uvotbadpix on $rawFile"); $ubadpix->params({infile => $rawFile, badpixlist => 'CALDB', outfile => $badFile, clobber => 'yes', chatter => 3}) ->run(); $log->entry("running uvotmodmap on $rawFile"); $umodmap->params({infile => $rawFile, badpixfile => $badFile, outfile => $mod8File, nsig => 3, clobber => 'yes', ncell => 2, # TODO: keeps run times for 2048^2 images to minutes # instead of hours chatter => 3}) ->run(); # uflatfield will go here $uflatfield->params({ infile => $mod8File, clobber => 'yes', infile => $rawFile, outfile => $corrFile}) ->run(); # unlink($mod8File); ######################################### # determine the DET coordinate filename ######################################### my $detFile = $filename->corresponding("rawimage", "detimage", $rawFile); $log->entry("converting $rawFile to $detFile"); ###################### # do the conversion ###################### $swiftxform->params({ infile => $mod8File, outfile => $detFile, }) ->run(); unlink($mod8File); } # end of loop over files } # end of image_raw_to_det method sub aspect_correct_sky_images { my $self = shift; my $log =$self->log(); my $filename=$self->filename(); my $procpar =$self->procpar(); my $jobpar =$self->jobpar(); my $jobtitle = $jobpar->read('job_title'); $log->entry('Performing aspect correction on UVOT images'); my $attfile = $filename->get('attitude', 's'); my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if(-f $uattfi){ $jumpatt = $uattfi; } $attfile = $jumpatt if -f $jumpatt; if (not $attfile) { $log->error(2, 'no attitude file available'); return; } my @skyFiles = $filename->get('filterimg', 'uvot', '*', '*'); if (not @skyFiles) { $log->entry("no sky images to correct"); return; } my $infile = 'uvotskycorr.files'; my $fh = FileHandle->new($infile, 'w'); if (not $fh) { $log->error(2, "unable to create $infile [$!]"); return; } foreach my $name (@skyFiles) { $fh->print($name . "\n"); } $fh->close; my $corrfile = $filename->get('hk', 'u', 'ac', 0); my $catspec = $procpar->read('starcatalog'); my $orgcor = undef; if(-e $corrfile and $jobtitle eq 'FINAL FOR ARCHIVE'){ $orgcor = $corrfile; $corrfile .= '_new'; } unlink($corrfile); $log->entry('finding corrections'); my $find = Util::HEAdas->new('uvotskycorr')->is_script(1); $find->params({ skyfile => '@' . $infile, what => 'ID', outfile => $corrfile, corrfile => 'NONE', attfile => $attfile, catspec => $catspec, starid => 'n.reference=50 n.observation=30 max.rate=1000', chatter => 5 }) ->run; unlink $infile; return if not -f $corrfile; $log->entry('applying corrections'); $fh = FileHandle->new($infile, 'w'); if (not $fh) { $log->error(2, "unable to create $infile [$!]"); return; } foreach my $name (@skyFiles) { $fh->print($name . "\n") if -f $name; $name =~ s/_sk([\._])/_ex$1/; $fh->print($name . "\n") if -f $name; } $fh->close; my $apply = Util::HEAdas->new('uvotskycorr')->is_script(1); $apply->params({ skyfile => '@' . $infile, what => 'SKY', outfile => 'NONE', corrfile => $corrfile, attfile => $attfile, catspec => $catspec, chatter => 5, }) ->run; ####################################################################### # Apply aspect corrections to attitude file. First check if there are # any corrections to apply. ####################################################################### my $corrfits = Util::FITSfile->new($corrfile)->cols('ASPCORR'); my @corrections = $corrfits->table(); if( grep(/1/, @corrections)){ $log->entry('applying corrections to attitude file'); my $scatt = $filename->get('attitude', 's'); my $jumpatt = $filename->get('attcorr', 'p'); my $uvotatt = $filename->get('attcorr', 'u'); # if(-f $uvotatt){ # rename $uvotatt, $uvotatt.'tmp'; # $attfile = $uvotatt.'tmp'; # } if (!-f $uvotatt) { my $att_corr = Util::HEAdas->new('uvotattcorr')->is_script(1); $att_corr->params({attfile => $attfile, corrfile => $corrfile, outfile => $uvotatt, clobber => 'yes'}) ->run(); } if (defined $orgcor and -e $orgcor) { my $fits = Util::FITSfile->new($orgcor); my $nhdus = $fits->nhdus(); my $Nname = 'ASPCORR'.$nhdus; Util::HEAdas->new('ftappend') ->params({ infile => "$corrfile\[ASPCORR\]", outfile => "$orgcor", }) ->run(); Util::HEAdas->new('fthedit') ->params({ infile => $orgcor.'['.$nhdus.']', keyword => 'EXTNAME', operation => 'add', value => $Nname }) ->run(); unlink $corrfile; } my $sattfits = Util::FITSfile->new($scatt, 0); my $attstatus = $sattfits->keyword('ATTSTATU'); # if( $att_corr->had_error() ){ # $attstatus += 2; # unlink $uvotatt; # rename $uvotatt.'tmp', $uvotatt; # }else{ # unlink $uvotatt.'tmp'; $attstatus += 1; my $flags = 111; $flags = 101 unless -f $jumpatt; if(-f $uvotatt){ my $uattfits = Util::FITSfile->new($uvotatt); for(my $i=0; $i<3; $i++){ $uattfits->ext($i); $uattfits->keyword('ATTFLAG', "'$flags'", 'Orgin of attitude information'); } } # } $sattfits->keyword('ATTSTATU', $attstatus, 'Status of corrected attitude files'); } } # end of aspect_correct_sky_images method ############################################################################### # generate exposure maps ############################################################################### sub create_exposure_maps { my $self=shift; my $log = $self->log(); my $filename= $self->filename(); my $procpar = $self->procpar(); my $jobpar = $self->jobpar(); $log->entry("Creating exposure maps"); ###################################################### # make sure there is an attitude file ###################################################### my $attitude = $filename->get('attitude', 's'); my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if (-f $uattfi) { $jumpatt = $uattfi; } $attitude = $jumpatt if -f $jumpatt; unless(-f $attitude) { $log->entry("No attitude data available - can't make exposure maps"); return; } ######################################### # create the swiftxform tool ######################################### ## attdelta => 100, my $uexpmap = Util::HEAdas->new('uvotexpmap') ->params({attfile => $attitude, teldeffile => 'CALDB', method => 'SHIFTADD', attdelta => 0.1, aberration => 'no', chatter => 4}); my $aspect_follow_file = $filename->get('hk', 'uvot', 'af', 0); if( -f $aspect_follow_file){ $uexpmap->params({trackfile => $aspect_follow_file}); }else{ $uexpmap->params({trackfile => 'NONE'}); } ###################### # loop over files ###################### foreach my $skyFile ($filename->get('filterimg', 'uvot', '*', '*') ) { # exclude grism files next if $filename->is_grism($skyFile, 'skyimage'); ######################################### # determine the exposure map filename ######################################### my ($expFile, $maskFile) = ($skyFile) x 2; $expFile =~ s/_sk([\._])/_ex$1/; $maskFile =~ s/_sk([\._])/_mk$1/; $log->entry("building exposure map $expFile for $skyFile"); my $badFile = $skyFile; $badFile =~ s/_sk([\._])/_bp$1/; if (not -f $badFile) { $log->entry("no bad pixel file for $skyFile"); $badFile = 'NONE'; } ####################### # make the exposure map ####################### $uexpmap->params({ infile => $skyFile, badpixfile => $badFile, outfile => $expFile, maskfile => $maskFile }) ->is_script(1) ->run(); } # end of loop over files } # end of create_exposure_maps method sub tag_round_images { my $self=shift; my $log = $self->log(); my $filename= $self->filename(); my $procpar = $self->procpar(); my $jobpar = $self->jobpar(); my $seq = $jobpar->read('sequence'); my ($code, $extname); if($seq%1000==992 || $seq%1000==0){ $log->entry("Find and tag images from the settling exposure."); my $catfile = $filename->get('hk', 'uvot', 'ct', 0); my $catfits = Util::FITSfile->new($catfile); my $crows = $catfits->nrows(); my %catvals; foreach my $col ('EXPID', 'MODEID', 'ESTART', 'ESTOP', 'FILTER'){ $catfits->cols($col); my @temp = $catfits->table(); $catvals{$col} = [@temp]; } my ($settle_id, $filter); for(my $row=0; $row<$crows; $row++){ my $cat_mode_id = $catvals{MODEID}->[$row]; next if $cat_mode_id eq 'NULL'; my $mode = $Subs::UvotNames::modeNames->[$cat_mode_id]; last if $mode eq 'Image/Event'; my $explen = $catvals{ESTOP}->[$row] - $catvals{ESTART}->[$row]; next unless ($mode eq 'Event' && $explen < 12); $settle_id = $catvals{EXPID}->[$row]; $filter = $catvals{FILTER}->[$row]; last; } unless($settle_id){ $log->entry("No settling exposure found."); }else{ $log->entry("Found settling exposure $settle_id, filter " . $Subs::UvotNames::filterNames->[$filter] .'.'); $code = $Subs::UvotNames::filterCodes->[$filter]; $extname = $code . $settle_id . 'E'; } } my @im_real_files = ($filename->get('filterimg', 'uvot', '', '*'), $filename->get('filterexp', 'uvot', '', '*')); my @imfiles = @im_real_files; push @imfiles, $filename->get('rawimage', 'uvot', '*', '*'); my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if (-f $uattfi) { $jumpatt = $uattfi; } foreach my $imfile (@imfiles){ my $imfits = Util::FITSfile->new($imfile); if( $code && $imfile=~/${code}_/ ){ my @hdus = $imfits->list_hdus(); if( grep(/${extname}/, @hdus) ){ $imfits->ext($extname); $imfits->keyword('OBS_MODE', 'SETTLING'); } } ######################## # Set the attitude flag ######################## if(-f $jumpatt){ $imfits->ext(0); if($jumpatt =~ /uat\.fits$/){ $imfits->keyword('ATTFLAG', "'111'", 'Orgin of attitude information'); } else { $imfits->keyword('ATTFLAG', "'110'", 'Orgin of attitude information'); } } } ############################################## # Round off real images so they compress well ############################################## my $tmp_round = 'round.tmp'; my $round = Util::HEAdas->new('swiftround') ->params({round => 0.01, outfile => $tmp_round}); foreach my $real_image (@imfiles){ $round->params({infile => $real_image}) ->run(); unless( $round->had_error() ){ unlink $real_image; rename $tmp_round, $real_image; } } } # end of tag_round_images sub sum_images { my $self=shift; my $log = $self->log(); my $filename= $self->filename(); my $procpar = $self->procpar(); my $jobpar = $self->jobpar(); $log->entry("Creating summed images"); ################################ # Create tool for summing images ################################ my $usumexp = Util::HEAdas->new('uvotimsum') ->params({clobber => 'yes'}); my $tmpFile = 'sum.tmp'; ######################################## # clean up any pre-existing summed files ######################################## foreach my $type (qw(skyimage expimage)) { foreach my $imageFile ($filename->get($type, 'uvot', '', '*')) { $log->entry("removing obsolete summed image $imageFile"); unlink($imageFile); } } ###################### # loop over files ###################### my %filter_to_sum = ('filterimg' => 'skyimage', 'filterexp' => 'expimage'); foreach my $type (keys %filter_to_sum) { $usumexp->params({exclude => 'DEFAULT'}); my @sumFiles; foreach my $filter_code (@$Subs::UvotNames::filterCodes){ next if ( $filter_code eq 'bl' || $filter_code eq 'un' ); my $filter = $filename->filter_name_for_code($filter_code); my @imageFiles = grep /u$filter_code/, $filename->get($type, 'uvot', '', '*'); next unless @imageFiles; my $sumFile = $filename->get($filter_to_sum{$type}, 'uvot', "${filter_code}sm", 0); push(@sumFiles, $sumFile); foreach my $imageFile (@imageFiles) { $log->entry("summing $imageFile to $sumFile"); my $maskFile = 'NONE'; unless( $filename->is_grism($imageFile, $type) ){ $maskFile = $imageFile; $maskFile =~ s/_[se][kx]([\._])/_mk$1/; } ###################### # do the summation ###################### $usumexp->params({ infile => $imageFile, outfile => $tmpFile, maskfile => $maskFile, expmap => ($type =~ /exp/ ? 'yes' : 'no'), }) ->is_script(1) ->run(); my $keys = "[col " . " #EXTNAME='$filter';" . " #FILTER='$filter';" . ']'; if (not -e $tmpFile) { $log->entry("no $filter summed image"); } elsif (not -e $sumFile) { Util::HEAdas->new('ftcopy') ->params({ infile => $tmpFile . "[col #EXTNAME='$filter']", outfile => $sumFile, copyall => 'yes', }) ->run(); }else { # append this result to the sum file Util::HEAdas->new('ftappend') ->params({ infile => "$tmpFile+1[col #EXTNAME='$filter']", outfile => $sumFile, }) ->run(); } } } $usumexp->params({exclude => 0}); my $sumFile = $filename->get($filter_to_sum{$type}, 'uvot', '', 0); if (-e $sumFile) { unlink $sumFile; } my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if (-f $uattfi) { $jumpatt = $uattfi; } foreach my $imageFile ($filename->get('filterimg', 'uvot', '*', '*')) { Util::HEAdas->new('ftcopy') ->params({ infile => $imageFile . "[0]", outfile => $sumFile, copyall => 'no', }) ->run(); Util::HEAdas->new('fthedit') ->params({ infile => $sumFile.'[0]', keyword => 'FILTER', operation => 'deleteall' }) ->run(); Util::HEAdas->new('fthedit') ->params({ infile => $sumFile.'[0]', keyword => 'HISTORY', operation => 'deleteall' }) ->run(); if (-f $jumpatt) { my $fits = Util::FITSfile->new($sumFile, 0); if($jumpatt =~ /uat\.fits$/){ $fits->keyword('ATTFLAG', "'111'", 'Orgin of attitude information'); } else { $fits->keyword('ATTFLAG', "'110'", 'Orgin of attitude information'); } } last; } foreach my $imageFile (@sumFiles) { next unless -f $imageFile; my (undef, $filter_code) = $filename->parse($imageFile, $filter_to_sum{$type}); $filter_code =~ s/sm$//; $log->entry("Filter code for $imageFile is $filter_code."); next if not $filter_code; ######################################################### # determine the summed image filename # by removing the filter code from sw\d{9}uff(sk|ex).img ######################################################### $log->entry("summing $imageFile to $sumFile"); my $filter = $filename->filter_name_for_code($filter_code); ###################### # do the summation ###################### $usumexp->params({ infile => $imageFile, outfile => $tmpFile, maskfile => 'NONE', expmap => ($type =~ /exp/ ? 'yes' : 'no'), }) ->is_script(1) ->run(); my $keys = "[col " . " #EXTNAME='$filter';" . " #FILTER='$filter';" . ']'; my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if(-f $uattfi){ $jumpatt = $uattfi; } if (not -e $tmpFile) { $log->entry("no $filter summed image"); }elsif (not -e $sumFile) { Util::HEAdas->new('ftcopy') ->params({ infile => $tmpFile . "[col #EXTNAME='$filter']", outfile => $sumFile, copyall => 'yes', }) ->run(); Util::HEAdas->new('fthedit') ->params({ infile => $sumFile.'[0]', keyword => 'EXTNAME', operation => 'delete' }) ->run(); if(-f $jumpatt){ my $fits = Util::FITSfile->new($sumFile, 0); if ($jumpatt =~ /uat\.fits$/) { $fits->keyword('ATTFLAG', "'111'", 'Orgin of attitude information'); } else { $fits->keyword('ATTFLAG', "'110'", 'Orgin of attitude information'); } } }else{ # append this result to the sum file Util::HEAdas->new('ftappend') ->params({ infile => "$tmpFile+1[col #EXTNAME='$filter']", outfile => $sumFile, }) ->run(); } unlink ($tmpFile, $imageFile); } # end of loop over files } # end of loop over types unlink $filename->get('maskimage', 'uvot', '*', '*'); } sub detect_sources { my $self=shift; my $log = $self->log(); my $filename= $self->filename(); my $procpar = $self->procpar(); my $jobpar = $self->jobpar(); $log->entry("Detecting sources"); my $jobtitle = $jobpar->read('job_title'); my $sexfile = $procpar->read('sexfile'); my $expopt = 'BETA'; if(!-e $sexfile){ $sexfile = 'DEFAULT'; } if($jobtitle ne 'FINAL FOR ARCHIVE'){ $sexfile = 'DEFAULT'; } if($sexfile ne 'DEFAULT'){ $expopt = 'ALPHA'; } ######################################## # clean up any pre-existing source lists ######################################## foreach my $catFile ($filename->get('srclist', 'uvot', '', '*')) { unlink($catFile); } ############################################################### # Create tools for detecting sources and determining magnitudes ############################################################### my $udetect = Util::HEAdas->new('uvotdetect'); my $tmpFile = 'sources.tmp'; my $catFile = undef; ############################# # loop over summed sky images ############################# foreach my $skyFile ($filename->get('skyimage', 'uvot', '', '*')) { $log->entry("detect: skyFile '$skyFile'"); my $expFile = $filename->corresponding('skyimage', 'expimage', $skyFile); if (!defined $expFile or !-e $expFile) { $log->entry("No exposure data for $skyFile - assuming constant exposure"); $expFile = 'NONE'; } if (not $catFile) { $catFile = $filename->corresponding('skyimage', 'srclist', $skyFile); } ######################################### # run uvotdetect on each extension ######################################### my $fitsFile = Util::FITSfile->new($skyFile, 0); my $nhdus = $fitsFile->nhdus(); my $jumpatt = $filename->get('attcorr', 'p'); my $uattfi = $filename->get('attcorr', 'u'); if(-f $uattfi){ $jumpatt = $uattfi; } for (my $hdu0 = 1; $hdu0 < $nhdus; ++$hdu0) { $fitsFile->ext($hdu0); my $filter = $fitsFile->keyword('EXTNAME'); next if $filter =~ /^g[uv]/; my $weightFile = $expFile; if ($expFile ne 'NONE') { $weightFile .= "[$filter]"; } if (-e $tmpFile) { unlink($tmpFile); } $udetect->params({ infile => "$skyFile+$hdu0", outfile => $tmpFile, # expfile => 'NONE', expfile => "$expFile+$hdu0", expopt => $expopt, sexfile => $sexfile, threshold => 3, # TODO: make parameter }) ->is_script(1) ->run(); if (-e $tmpFile) { if (not -e $catFile) { Util::HEAdas->new('ftcopy') ->params({ infile => $tmpFile . "[col #EXTNAME='$filter']", outfile => $catFile, copyall => 'yes', }) ->run(); if (-f $jumpatt) { my $fits = Util::FITSfile->new($catFile, 0); if ($jumpatt =~ /uat\.fits$/) { $fits->keyword('ATTFLAG', "'111'", 'Orgin of attitude information'); } else { $fits->keyword('ATTFLAG', "'110'", 'Orgin of attitude information'); } } } else { # append this result to the sum file Util::HEAdas->new('ftappend') ->params({ infile => "$tmpFile+1[col #EXTNAME='$filter']", outfile => $catFile, }) ->run(); } } # end checking if temp file exists } # end of loop over extensions } # end of loop over files unlink($tmpFile) if -e $tmpFile; } # end of detect_sources method ################################################################################ # plot images ################################################################################ sub plot_images { my $self = shift; my $log =$self->log(); my $filename=$self->filename(); my $procpar =$self->procpar(); my $jobpar =$self->jobpar(); my $max_image_dimen=1024; $log->entry("Plotting all images"); my $ximage = Util::Xanadu->new("ximage"); my $type; my $rate_img = "rate_image.tmp"; foreach $type ("detimage", "skyimage") { my $plot_type = $type; $plot_type =~ s/image/plot/; ###################################### # loop over summed images of this type ###################################### my $file; my $title = uc($jobpar->read('object')) .' SWIFT '. $jobpar->read('sequence') .' '; foreach my $file ($filename->get($type, 'uvot', '', '*')) { my $title2; my $fitsFile = Util::FITSfile->new($file, 0); my $nhdus = $fitsFile->nhdus(); my $grid=''; my $arith; if($type =~ /^sky/ ) { $grid = 'grid'; $arith = Util::HEAdas->new('farith') ->params({overflow => 'yes', null => 'yes', ops => 'DIV', outfil => $rate_img, clobber => 'yes'}); } my $date = $fitsFile->keyword('DATE-OBS'); $date =~ s/T.*$//; my $title1 = $title.$date; my $inst = $fitsFile->keyword('INSTRUME'); my $stat = Util::HEAdas->new('ftstat'); # ->params({centroid => 'no'}); ############################################## # for each extension (filter) create a graphic ############################################## for(my $hdu=1; $hdu<$nhdus; $hdu++) { my $inplot = "$file\[$hdu\]"; ################################################## # If a sky image, make a rate image and plot that ################################################## if($type =~ /^sky/ ){ my $infits = Util::FITSfile->new($file, $hdu); my $expmap = $file .'['. $infits->keyword('EXTNAME') .']'; # my $expmap = $file .'['. $hdu .']'; $expmap =~ s/_sk/_ex/; $arith->params({infil1 => $inplot, infil2 => $expmap}) ->run(); $inplot = $rate_img; } $fitsFile->ext($hdu); my $fitsFile2 = Util::FITSfile->new($inplot, 0); $fitsFile2->ext(0); my $naxis1 = $fitsFile2->keyword('NAXIS1'); my $naxis2 = $fitsFile2->keyword('NAXIS2'); my $xulm = int($naxis1/3.0); my $yulm = int($naxis2/3.0); my $xllm = int(2.0*$naxis1/3.0); my $yllm = int(2.0*$naxis2/3.0); my $strL = "\[$xulm:$yulm,$xllm:$yllm\]"; my $partfile = $inplot."\[0\]$strL"; $stat->params({infile => $partfile})->run(); my $stdo = $stat->stdout(); my $par = $stat->parfile(); my $mean = $par->read('mean'); $stat->params({infile => $inplot})->run(); my $par2 = $stat->parfile(); my $tmin = $par2->read('min'); my $tmax = $par2->read('max'); my $median = $self->getMedian($stdo); my $min = $median*3.0; my $max = $median*400; my @commands; push @commands, ('cpd /vgif', 'cey 2000'); push @commands, ('color/setcolor=0/namecolor=white'); push @commands, ('color/setcolor=1/namecolor=black'); my $filter = $fitsFile->keyword('EXTNAME'); my $exp = $fitsFile->keyword('EXPOSURE'); my $filter_code = $filename->filter_code_for_name($filter); $title2 = $inst .' '. $filter .' Exposure '. sprintf('%10d', int($exp)) .'s'; my $plotFile = $filename->get($plot_type, 'uvot', $filter_code, 0); ######################################## # determine if we need to bin the image ######################################## my $dimen = $fitsFile->keyword('NAXIS1'); my $bin=1; for($bin=1; $max_image_dimen*$bin<$dimen; $bin++) {} $log->entry("Plotting ${file}\[$hdu\] in $plotFile binned by $bin"); push @commands, ("read/fits/rebin=$bin \{$inplot\}", 'levels/num=64', 'cct/set invgray', "title \"$title1\"", "title/lower \"$title2\"", "disp/log/minlevel=". $min ."/maxlevel=". $max, $grid, 'scale', 'show'); push @commands, 'exit'; $ximage->script(@commands)->run(); rename 'pgplot.gif', $plotFile; } } } unlink $rate_img; } # end of plot_images method sub getMedian { my $self = shift; my ($str) = @_; my @tar = split/\n/, $str; my $median = undef; foreach my $l (@tar){ next if($l =~ /^\s*$/ or $l =~ /^\s*\=/); $l =~ s/^\s+//; if($l =~ /median\:/){ $median = (split/\:\s+/, $l)[-1]; last; } } return $median; } 1;