7 June 2005 Reprocessing our Chandra HETGS observation of DoAr 21 based on Marc's suggestion Have just installed CIAO v3.2.1 and CALDB v3.0.3 Working in: /home/cohen/doar21/200231/secondary/ tarred up all the files that were originally in there Following the "Obtain Grating Spectra from HETG/ACIS-S Data" thread(s) at: http://cxc.harvard.edu/ciao/threads/spectra_hetgacis/ - gunzipped the evt1 file: acisf03761_000N001_evt1.fits.gz - checked CALDB version used to process this file: [cohen@ruby secondary]$ dmkeypar acisf03761_000N001_evt1.fits CALDBVER echo+ Warning: found possibly out of date user parameter file, renamed to /home/cohen/cxcds_param/dmkeypar_20050607.15:29:27.par 2.21 - and the processing version: [cohen@ruby secondary]$ dmkeypar acisf03761_000N001_evt1.fits ASCDSVER echo+ 6.13.0 OK, this is lower than 7.4.0, and the data perparation information says: If your data was processed with a software version lower than DS 7.4.0, it is recommended that you remove the acis_detect_afterglow correction and reprocess the data with a new bad pixel file before running this thread. The thread for doing this is at: http://cxc.harvard.edu/ciao/threads/acisdetectafterglow/ Note that the thread says: "After running this thread, users should complete the Identify ACIS Hot Pixels and Cosmic Ray Afterglows thread." But wait, notes linked from the thread say: "Acis_detect_afterglow rejects 3-5% of valid source photons from diffracted spectra. The rejection is not "gray". For line sources, the rejections occur in strong lines. It may not be gray for continuum spectra, either: John's pileup tests on ObsID 105 had very localized residual features." and "Do NOT run acis_detect_afterglow for any grating (HETG/ACIS-S or LETG/ACIS-S) data. Acis_detect_afterglow should be removed from the grating pipelines." OK, so back to the original data preparation thread: Generate A New Level=1.5 Event File 1. Get position of zero-order image (tgdetect) To find the zero-order location, the tool tgdetect is run: [cohen@ruby secondary]$ punlearn tgdetect [cohen@ruby secondary]$ pset tgdetect infile=acisf03761_000N001_evt1.fits [cohen@ruby secondary]$ pset tgdetect outfile=acisf03761_000N001_src1a.fits [cohen@ruby secondary]$ tgdetect OK...(with spurious error message "double" - see bug page) Overlayed the source list onto the image, using DS9, it looks right (to the naked eye, the source circle does NOT look perfectly centered on the zeroth order spectrum). Note: I saved a screen shot into doar21_reprocessing.ppt Next: 2. Get region mask (tg_create_mask) 3. Run tg_resolve_events *I'm a little confused about the various input files. I copied in the one asol1 file from the "primary" directory. Note: I set the eventdef parameter to ")stdlev1_ACIS" as the thread indicates. The program chugs for a few seconds, then gives the unilluminating error message that it couldn't create the output file. 8 July 05 - resuming this project... We think the problem might have been caused by my not restarting my vnc session after all the software upgrades, etc. Marc reported that he ran this complete process just fine on the very same files on his machine. Oh my God! The problem was that a file existed, but with the .gz on it, and ciao wouldn't overwrite that. So, new file is named: acisf03761_000N002_evt1a.fits Now, Generate a new level 2 event file: 1. Apply grade/status filters (dmcopy) [cohen@ruby secondary]$ punlearn dmcopy [cohen@ruby secondary]$ dmcopy "acisf03761_000N002_evt1a.fits[EVENTS][grade=0,2,3,4,6,status=0]" acisf03761_000N001_flt1_evt1a.fits opt=all 2. Apply GTI filters (dmcopy) [cohen@ruby secondary]$ punlearn dmcopy [cohen@ruby secondary]$ dmcopy "acisf03761_000N001_flt1_evt1a.fits[EVENTS][@acisf03761_000N001_flt1.fits][cols -phas]" acisf03761_000N001_evt2.fits opt=all Run destreak [cohen@ruby secondary]$ punlearn destreak [cohen@ruby secondary]$ pset destreak infile=acisf03761_000N001_evt2.fits [cohen@ruby secondary]$ pset destreak outfile=acisf03761_000N001_dstrk_evt2.fits *check ppt file for screenshot of ds9 showing destreaked evt2 file image Extract a Grating Spectrum (tgextract) [cohen@ruby secondary]$ punlearn tgextract [cohen@ruby secondary]$ pset tgextract infile=acisf03761_000N001_dstrk_evt2.fits [cohen@ruby secondary]$ pset tgextract outfile=acisf03761_000N001_dstrk_pha2.fits [cohen@ruby secondary]$ tgextract Input event file (output event file from L1.5 processing) (acisf03761_000N001_dstrk_evt2.fits): If typeII, enter full output file name or '.'; if typeI, enter output rootname (acisf03761_000N001_dstrk_pha2.fits): Input ancillary response file name (none): Input redistribution file name (none): Source ID's to process: 'all', comma list, @file (all): Grating parts to process: HETG, HEG, MEG, LETG, header_value (HETG|HEG|MEG|LETG|header_value) (header_value): Grating diffraction orders to process: 'default', comma list, range list, @file (default): Output file type: typeI (single spectrum) or typeII (multiple spectra) (pha_typeI|pha_typeII) (pha_typeII): Next, make grating rmfs: Create Grating RMFs for ACIS-S Observations http://cxc.harvard.edu/ciao/threads/mkgrmf_aciss/ [cohen@ruby secondary]$ punlearn mkgrmf [cohen@ruby secondary]$ pset mkgrmf order=1 [cohen@ruby secondary]$ pset mkgrmf grating_arm=HEG [cohen@ruby secondary]$ pset mkgrmf outfile=doar21_heg_p1.rmf [cohen@ruby secondary]$ pset mkgrmf obsfile="acisf03761_000N001_dstrk_pha2.fits[SPECTRUM]" [cohen@ruby secondary]$ pset mkgrmf regionfile=acisf03761_000N001_dstrk_pha2.fits [cohen@ruby secondary]$ pset mkgrmf detsubsys=ACIS-S3 [cohen@ruby secondary]$ pset mkgrmf wvgrid_arf=compute [cohen@ruby secondary]$ pset mkgrmf wvgrid_chan=compute [cohen@ruby secondary]$ pset mkgrmf clobber=no [cohen@ruby secondary]$ mkgrmf NOTE: this is for the +1 HEG grmf... repeat for all orders and both gratings OK... first we need to make sure that ardlib.par points at the right bad pixel file...and apparently, we also need to remake the bad pixel file. (then redo the above mkgrmf) - Marc pointed out that the rmf's (as opposed to grmfs are just copied from caldb directories, and also that they now have separate ones for -1 and +1). Identify ACIS Hot Pixels and Cosmic Ray Afterglows http://cxc.harvard.edu/ciao/threads/acishotpixels/ *need to install the 3.2.2 patch 28July05 - OK, patch installed. Here's the sequence of things we need to do: http://cxc.harvard.edu/ciao/threads/mkgrmf_aciss/ !create grating rmfs; but first: http://cxc.harvard.edu/ciao/threads/badpix/ !make sure ardlib.par points at the right bad pixels file; but first: http://cxc.harvard.edu/ciao/threads/acishotpixels/ !make a new bad pixels file; but first: http://cxc.harvard.edu/ciao/threads/acisdetectafterglow/ !remove the old (standard data processing) afterglow correction: OK, doing the remove afterglow - I'm a little concerned because the input is evt1, and we've already used this (and some other things) to make a new evt2 file. But going ahead... ------ ! ! ! ------ 3 Aug 05: in new folder secondary/reproc/ - copied in evt1.fits file - used dmkeypar to find out that processing was with version 6.13.0 - removing afterglow correction: http://cxc.harvard.edu/ciao/threads/acisdetectafterglow/ creating copy that shows the bits flagged as afterglows: [cohen@ruby reproc]$ dmcopy "acisf03761_000N001_evt1.fits[exclude status=xxxxxxxxxxxx0000xxxxxxxxxxxxxxxx]" acis_afterglows_evt1.fits the source (zeroth order, anyway) is clearly visible. Indicating that source counts have been flagged as afterglow events (incorrectly). [cohen@ruby reproc]$ punlearn dmtcalc [cohen@ruby reproc]$ pset dmtcalc infile=acisf03761_000N001_evt1.fits [cohen@ruby reproc]$ pset dmtcalc outfile=acisf03761_000N001_reset_evt1.fits [cohen@ruby reproc]$ pset dmtcalc expression="status=status,status=X15F,status=X14F,status=X13F,status=X12F" [cohen@ruby reproc]$ dmtcalc *Note the name of the new file (with "reset" put in the middle); this is now our working evt1 file. Took a long time to run, but eventually completed without any error messages. acisf03761_000N001_reset_evt1.fits is our new working evt1 file (it's almost identical in size to the original one) At the end of this thread it says: At this point, users should procede to the Identify ACIS Hot Pixels and Cosmic Ray Afterglows thread to create a new bad pixel file for their data which contains more accurate afterglow information. After that, the data is reprocessed with acis_process_events and a new level=2 file must be created. For an overview of the full reprocessing thread, see the ACIS Data Preparation Analysis Guide. So, moving on - http://cxc.harvard.edu/ciao/threads/acishotpixels/ Copied the (old) bpix1.fits file in from the primary/ directory, and copied in the msk1.fits file (after gunzipping it) from the secondary/ directory. Retrieving the bias and block parameter files from the web, using our OBSID=3761 The name of your tar file is: B8249.tar You may retrieve it from our ftp site cxc.harvard.edu, from the directory /pub/arcftp/webstage/B8249 in the next eight hours. These actually constitute 7 files - six bias files and one pbk0 file [cohen@ruby reproc]$ punlearn acis_run_hotpix [cohen@ruby reproc]$ pset acis_run_hotpix infile=acisf03761_000N001_reset_evt1.fits [cohen@ruby reproc]$ pset acis_run_hotpix outfile=acisf03761_000N001_new_bpix1.fits [cohen@ruby reproc]$ pset acis_run_hotpix badpixfile=acisf03761_000N001_bpix1.fits [cohen@ruby reproc]$ pset acis_run_hotpix biasfile=@bias_files.lis [cohen@ruby reproc]$ pset acis_run_hotpix maskfile=acisf03761_000N001_msk1.fits [cohen@ruby reproc]$ pset acis_run_hotpix pbkfile=acisf168496622N001_pbk0.fits [cohen@ruby reproc]$ emacs bias_files.lis & Note: make the bias_files.lis - all six of the bias files So, new badpix file has been created: acisf03761_000N001_new_bpix1.fits Next, we apply this correction: [cohen@ruby reproc]$ pset acis_process_events badpixfile=acisf03761_000N001_new_bpix1.fits and run the acis_process_events task First, make sure ardlib has right bpix file: http://cxc.harvard.edu/ciao/threads/badpix/index.html#punlearn [cohen@ruby reproc]$ punlearn ardlib [cohen@ruby reproc]$ acis_set_ardlib acisf03761_000N001_new_bpix1.fits Updated ardlib parameter file: /home/cohen/pfiles/ardlib.par AXAF_ACIS0_BADPIX_FILE -> CALDB AXAF_ACIS1_BADPIX_FILE -> CALDB AXAF_ACIS2_BADPIX_FILE -> CALDB AXAF_ACIS3_BADPIX_FILE -> CALDB AXAF_ACIS4_BADPIX_FILE -> /home/cohen/doar21/200231/secondary/reproc/acisf03761_000N001_new_bpix1.fits[BADPIX4] AXAF_ACIS5_BADPIX_FILE -> /home/cohen/doar21/200231/secondary/reproc/acisf03761_000N001_new_bpix1.fits[BADPIX5] AXAF_ACIS6_BADPIX_FILE -> /home/cohen/doar21/200231/secondary/reproc/acisf03761_000N001_new_bpix1.fits[BADPIX6] AXAF_ACIS7_BADPIX_FILE -> /home/cohen/doar21/200231/secondary/reproc/acisf03761_000N001_new_bpix1.fits[BADPIX7] AXAF_ACIS8_BADPIX_FILE -> /home/cohen/doar21/200231/secondary/reproc/acisf03761_000N001_new_bpix1.fits[BADPIX8] AXAF_ACIS9_BADPIX_FILE -> /home/cohen/doar21/200231/secondary/reproc/acisf03761_000N001_new_bpix1.fits[BADPIX9] OK, so now we go back to the ACIS-S HETGS data preparation thread: http://cxc.harvard.edu/ciao/threads/spectra_hetgacis/ We'll be creating a new evt1.5 file and then a evt2 file copying in the two additional files: [cohen@ruby secondary]$ cp pcadf168496624N001_asol1.fits reproc/. [cohen@ruby secondary]$ cp acisf03761_000N001_flt1.fits reproc/. Generate A New Level=1.5 Event File 1. Get position of zero-order image (tgdetect) To find the zero-order location, the tool tgdetect is run: *problem copying temporary tgdetect.par from current working dir [cohen@ruby reproc]$ setenv PFILES ${PFILES}:. [cohen@ruby reproc]$ tgdetect Input L1 event file (acisf03761_000N001_reset_evt1.fits): Input source position(s) file from previous OBI or NONE (NONE): Output source position(s) file name (acisf03761_000N001_reset_src1a.fits): # DMCOPY (CIAO3.2): Filter data type mismatch for DOUBLE Note: there's a bug report about this message (don't worry about it) Note also, that the tgdetect.par file (and temp files) are in the working directory, not in cxcds_param/ Viewing the evt1 file with the source overlaid: 2. Get region mask (tg_create_mask) The location of the HEG and MEG "arms" needs to be found next, via the tool tg_create_mask. This tool creates a region file, acis_evt1_L1a.fits, that will be used later to mask the image: [cohen@ruby reproc]$ punlearn tg_create_mask [cohen@ruby reproc]$ pset tg_create_mask infile=acisf03761_000N001_reset_evt1.fits [cohen@ruby reproc]$ pset tg_create_mask outfile=acisf03761_000N001_reset_evt1_L1a.fits [cohen@ruby reproc]$ pset tg_create_mask input_pos_tab=acisf03761_000N001_reset_src1a.fits [cohen@ruby reproc]$ tg_create_mask 3. Run tg_resolve_events The tool tg_resolve_events is now used to assign grating events to spectral orders, using the detector energy resolution for order separation: [cohen@ruby reproc]$ punlearn tg_resolve_events [cohen@ruby reproc]$ pset tg_resolve_events infile=acisf03761_000N001_reset_evt1.fits [cohen@ruby reproc]$ pset tg_resolve_events outfile=acisf03761_000N001_reset_evt1a.fits [cohen@ruby reproc]$ pset tg_resolve_events regionfile=acisf03761_000N001_reset_evt1_L1a.fits [cohen@ruby reproc]$ pset tg_resolve_events acaofffile=pcadf168496624N001_asol1.fits [cohen@ruby reproc]$ pset tg_resolve_events eventdef=")stdlev1_ACIS" [cohen@ruby reproc]$ tg_resolve_events Generate A New Level=2 Event File 1. Apply grade/status filters (dmcopy) Filter for bad grades (using ASCA grades) and for bit status: [cohen@ruby reproc]$ punlearn dmcopy [cohen@ruby reproc]$ dmcopy "acisf03761_000N001_reset_evt1a.fits[EVENTS][grade=0,2,3,4,6,status=0]" acisf03761_000N001_reset_flt1_evt1a.fits opt=all 2. Apply GTI filters (dmcopy) The Good Time Intervals (GTIs) supplied by the pipeline now need to be applied. We simultaneously eliminate an unnecessary column from the output: [cohen@ruby reproc]$ punlearn dmcopy [cohen@ruby reproc]$ dmcopy "acisf03761_000N001_reset_flt1_evt1a.fits[EVENTS][@acisf03761_000N001_flt1.fits][cols -phas]" acisf03761_000N001_reset_flt1_evt2.fits opt=all Run destreak There is a flaw in the serial readout of the ACIS chips, causing a significant amount of charge to be randomly deposited along pixel rows as they are read out. Although not much coincidence is expected for low-rate/low-exposure sources, ACIS-S4 (ccd_id=8) is significantly affected by this problem. The destreak tool detects coincidence of events in adjacent pixels along a row, flags probable streak events, and removes them. For details on how the tool works, see the Destreaking ACIS Data why topic. The streak events are not removed by standard grade, status, or bad pixel filtering, and are therefore still present in the level=2 event file. this image [Link to Image 2] compares the original level=2 event file (left) and the destreaked version (right). The ccd_id parameter should always be set to 8 (the default) to ensure that the destreak algorithm will be applied to the ACIS-S4 chip only. The algorithm should not be applied to chips that do not show streaks like those on ACIS-S4; see ahelp destreak for more information. [cohen@ruby reproc]$ punlearn destreak [cohen@ruby reproc]$ pset destreak infile=acisf03761_000N001_reset_flt1_evt2.fits [cohen@ruby reproc]$ pset destreak outfile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits [cohen@ruby reproc]$ destreak Extract a Grating Spectrum (tgextract) The CIAO tool tgextract produces a PHA2 spectrum file from the level=2 event file: [cohen@ruby reproc]$ punlearn tgextract [cohen@ruby reproc]$ pset tgextract infile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits [cohen@ruby reproc]$ pset tgextract outfile=acisf03761_000N001_reset_flt1_dstrk_pha2.fits [cohen@ruby reproc]$ tgextract Input event file (output event file from L1.5 processing) (acisf03761_000N001_reset_flt1_dstrk_evt2.fits): If typeII, enter full output file name or '.'; if typeI, enter output rootname (acisf03761_000N001_reset_flt1_dstrk_pha2.fits): Input ancillary response file name (none): Input redistribution file name (none): Source ID's to process: 'all', comma list, @file (all): Grating parts to process: HETG, HEG, MEG, LETG, header_value (HETG|HEG|MEG|LETG|header_value) (header_value): Grating diffraction orders to process: 'default', comma list, range list, @file (default): Output file type: typeI (single spectrum) or typeII (multiple spectra) (pha_typeI|pha_typeII) (pha_typeII): NOTE: names: I've just appended name fragments (e.g. 'dstrk') each time a process is performed on a file. I will make symbolic links to the "final" or "good" evt2 and pha2 files, so it's clear which are the ones to use. We could also make a hard copy for the purpose of performing analysis... Summary This thread is complete; the PHA2 grating spectrum file is named acis_pha2.fits. You should now proceed to the Create Grating RMFs for ACIS-S Observations thread. Run mkgrmf Here we create a +1 order HEG gRMF with standard grids: Aside: About the calibration files The calibration files used by mkgrmf (LSFPARM files) for the HETG grating arm (HEG, MEG) are labeled as +1 and -1: unix% ls -1 $CALDB/data/chandra/tel/grating/hetg/cpf/lsf/ acisheg-1D1999-07-22lsfparmN0003.fits acisheg-1D1999-07-22lsfparmN0004.fits acisheg1D1999-07-22lsfparmN0003.fits acisheg1D1999-07-22lsfparmN0004.fits acismeg-1D1999-07-22lsfparmN0003.fits acismeg-1D1999-07-22lsfparmN0004.fits acismeg1D1999-07-22lsfparmN0003.fits acismeg1D1999-07-22lsfparmN0004.fits They scale nicely at higher orders, however, allowing one to build a fair representation of 2nd or 3rd order gRMFs based on 1st order LSF model. [cohen@ruby reproc]$ punlearn mkgrmf [cohen@ruby reproc]$ pset mkgrmf order=1 [cohen@ruby reproc]$ pset mkgrmf grating_arm=HEG [cohen@ruby reproc]$ pset mkgrmf outfile=heg_p1.rmf [cohen@ruby reproc]$ pset mkgrmf obsfile="acisf03761_000N001_reset_flt1_dstrk_pha2.fits[SPECTRUM]" [cohen@ruby reproc]$ pset mkgrmf regionfile=acisf03761_000N001_reset_flt1_dstrk_pha2.fits [cohen@ruby reproc]$ pset mkgrmf detsubsys=ACIS-S3 [cohen@ruby reproc]$ pset mkgrmf wvgrid_arf=compute [cohen@ruby reproc]$ pset mkgrmf wvgrid_chan=compute [cohen@ruby reproc]$ pset mkgrmf clobber=no [cohen@ruby reproc]$ mkgrmf Output File Name (heg_p1.rmf): Enter ARF side wavelegth grid [angstroms] (compute): Enter channel-side wavelegth grid [angstroms] (compute): Enter Grating order (1): Name of fits file with obs info (acisf03761_000N001_reset_flt1_dstrk_pha2.fits[SPECTRUM]): File containing extraction region (acisf03761_000N001_reset_flt1_dstrk_pha2.fits): SrcID (1): Enter RMF threshold (1e-06): Verbosity (0:5) (0): Detector Name (e.g., ACIS-S3) (ACIS-S3): Enter Grating Arm (HEG|MEG|LEG|NONE) (HEG): Then - Run the tool for every order and arm you wish to model, changing the order and grating_arm parameters accordingly. So, to do the minus 1st order grmf: [cohen@ruby reproc]$ mkgrmf Output File Name (heg_p1.rmf): heg_m1.rmf Enter ARF side wavelegth grid [angstroms] (compute): Enter channel-side wavelegth grid [angstroms] (compute): Enter Grating order (1): -1 Name of fits file with obs info (acisf03761_000N001_reset_flt1_dstrk_pha2.fits[SPECTRUM]): File containing extraction region (acisf03761_000N001_reset_flt1_dstrk_pha2.fits): SrcID (1): Enter RMF threshold (1e-06): Verbosity (0:5) (0): Detector Name (e.g., ACIS-S3) (ACIS-S3): Enter Grating Arm (HEG|MEG|LEG|NONE) (HEG): Now, making the -1 and +1 MEG rmfs: [cohen@ruby reproc]$ mkgrmf Output File Name (heg_m1.rmf): meg_m1.rmf Enter ARF side wavelegth grid [angstroms] (compute): Enter channel-side wavelegth grid [angstroms] (compute): Enter Grating order (-1): Name of fits file with obs info (acisf03761_000N001_reset_flt1_dstrk_pha2.fits[SPECTRUM]): File containing extraction region (acisf03761_000N001_reset_flt1_dstrk_pha2.fits): SrcID (1): Enter RMF threshold (1e-06): Verbosity (0:5) (0): Detector Name (e.g., ACIS-S3) (ACIS-S3): Enter Grating Arm (HEG|MEG|LEG|NONE) (HEG): MEG [cohen@ruby reproc]$ mkgrmf Output File Name (meg_m1.rmf): meg_p1.rmf Enter ARF side wavelegth grid [angstroms] (compute): Enter channel-side wavelegth grid [angstroms] (compute): Enter Grating order (-1): 1 Name of fits file with obs info (acisf03761_000N001_reset_flt1_dstrk_pha2.fits[SPECTRUM]): File containing extraction region (acisf03761_000N001_reset_flt1_dstrk_pha2.fits): SrcID (1): Enter RMF threshold (1e-06): Verbosity (0:5) (0): Detector Name (e.g., ACIS-S3) (ACIS-S3): Enter Grating Arm (HEG|MEG|LEG|NONE) (MEG): OK, the last thing to do before actually analyzing the data is to make grating ARFs: http://asc.harvard.edu/ciao/threads/mkgarf_hetgacis/ [cohen@ruby reproc]$ grep version_ 'which fullgarf' grep: which fullgarf: No such file or directory Looks like I don't even have the script... Following instructions for installing the scripts at: http://asc.harvard.edu/ciao/download/scripts/README_CIAO_scripts The scripts themselves are available as a single tarball at: http://asc.harvard.edu/ciao/download/scripts/#scripts OK, despite the message above, fullgarf is in the scripts directory: /opt/ciao/contrib/bin and it was installed on May 6, 05; the website says it hasn't been updated since 2003, so I am NOT going to install the scripts tarball. But it is sitting in /home/cohen/chandra/ if I want to use it later. Note: ardlib.par seems to be pointing at the right contamination files. See: http://asc.harvard.edu/ciao/threads/aciscontam/index.html#setparam [cohen@ruby reproc]$ plist ardlib |grep CONTAM AXAF_ACIS0_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS1_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS2_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS3_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS4_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS5_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS6_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS7_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS8_CONTAM_FILE = CALDB Enter ACIS Contamination File AXAF_ACIS9_CONTAM_FILE = CALDB Enter ACIS Contamination File Run fullgarf To make sure that the gARF file is correct, we need to set the maskfile parameter for mkgarf. Currently, this cannot be done via fullgarf; it must be specified in the mkgarf parameter file before running the script: unix% pset mkgarf maskfile=acisf00459_000N002_msk1.fits Here's my issue of this command: [cohen@ruby reproc]$ pset mkgarf maskfile=acisf03761_000N001_msk1.fits Now we have all the information needed to run fullgarf. We will have to run the script twelve times, once for each row in the PHA file. But we're only doing first order garfs right now (so 4 times); that's all the rmfs we have. [cohen@ruby reproc]$ pset fullgarf phafile=doar21_pha2.fits [cohen@ruby reproc]$ pset fullgarf pharow=1 [cohen@ruby reproc]$ pset fullgarf evtfile=doar21_evt2.fits [cohen@ruby reproc]$ pset fullgarf asol=pcadf168496624N001_asol1.fits [cohen@ruby reproc]$ pset fullgarf pharow=3 [cohen@ruby reproc]$ pset fullgarf engrid="grid(heg_m1.rmf[cols ENERG_LO,ENERG_HI])" [cohen@ruby reproc]$ pset fullgarf dtffile=")evtfile" [cohen@ruby reproc]$ pset fullgarf badpix=acisf03761_000N001_new_bpix1.fits [cohen@ruby reproc]$ pset fullgarf rootname=acisf03761 [cohen@ruby reproc]$ fullgarf Error in parameter file: line 5: field Prompt Could not find the parameter file! (Looked for ) Aborting. Hmmm...it's like CIAO can't find the script (or the param file for it)... OK, reentering params and trying again: [cohen@ruby reproc]$ punlearn fullgarf [cohen@ruby reproc]$ pset fullgarf pharow=3 [cohen@ruby reproc]$ fullgarf Will use /home/cohen/cxcds_param/fullgarf.par for the parameter file. /home/cohen/cxcds_param/fullgarf.par contains 11 parameters . . . Input PHA file (Type I or II) (): acisf03761_000N001_reset_flt1_dstrk_pha2.fits Row in Type II PHA file (ignored if Type I) (3): Event file (): acisf03761_000N001_reset_flt1_dstrk_evt2.fits Aspect offsets file (): pcadf168496624N001_asol1.fits Energy grid spec (): "grid(heg_m1.rmf[cols ENERG_LO,ENERG_HI])" Dead time correction factor; ACIS->evt file; HRC -> dtf file ()evtfile): acisf03761_000N001_reset_flt1_dstrk_evt2.fits Bad pixel file; (filename|NONE|CALDB) (NONE): acisf03761_000N001_new_bpix1.fits Output rootname (): acisf03761 Getting the pha file type . . . Grating arm is HEG, order=-1 Source location is X=4131.1704101562, Y=4040.8210449219 Keyword existence (yes): Detector is ACIS Will run asphist for ccd_id= 4 5 6 7 asphist infile=pcadf168496624N001_asol1.fits[@acisf03761_000N001_reset_flt1_dstrk_evt2.fits[ccd_id=4]] outfile=acisf03761_ah4.fits evtfile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits dtffile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits asphist infile=pcadf168496624N001_asol1.fits[@acisf03761_000N001_reset_flt1_dstrk_evt2.fits[ccd_id=5]] outfile=acisf03761_ah5.fits evtfile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits dtffile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits asphist infile=pcadf168496624N001_asol1.fits[@acisf03761_000N001_reset_flt1_dstrk_evt2.fits[ccd_id=6]] outfile=acisf03761_ah6.fits evtfile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits dtffile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits asphist infile=pcadf168496624N001_asol1.fits[@acisf03761_000N001_reset_flt1_dstrk_evt2.fits[ccd_id=7]] outfile=acisf03761_ah7.fits evtfile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits dtffile=acisf03761_000N001_reset_flt1_dstrk_evt2.fits Finished processing aspect histograms for ccd_id= 4 5 6 7 Will run mkgarf for the same ccd_id list mkgarf detsubsys=ACIS-S0 order=-1 grating_arm=HEG outfile=acisf03761_S0_HEG_-1.fits asphistfile=acisf03761_ah4.fits[ASPHIST] engrid=grid(heg_m1.rmf[cols ENERG_LO osipfile=CALDB mode=hl verb=0 **** grid spec "grid(heg_m1.rmf[cols ENERG_LO" specifies an invalid grid mkgarf detsubsys=ACIS-S1 order=-1 grating_arm=HEG outfile=acisf03761_S1_HEG_-1.fits asphistfile=acisf03761_ah5.fits[ASPHIST] engrid=grid(heg_m1.rmf[cols ENERG_LO osipfile=CALDB mode=hl verb=0 **** grid spec "grid(heg_m1.rmf[cols ENERG_LO" specifies an invalid grid mkgarf detsubsys=ACIS-S2 order=-1 grating_arm=HEG outfile=acisf03761_S2_HEG_-1.fits asphistfile=acisf03761_ah6.fits[ASPHIST] engrid=grid(heg_m1.rmf[cols ENERG_LO osipfile=CALDB mode=hl verb=0 **** grid spec "grid(heg_m1.rmf[cols ENERG_LO" specifies an invalid grid mkgarf detsubsys=ACIS-S3 order=-1 grating_arm=HEG outfile=acisf03761_S3_HEG_-1.fits asphistfile=acisf03761_ah7.fits[ASPHIST] engrid=grid(heg_m1.rmf[cols ENERG_LO osipfile=CALDB mode=hl verb=0 **** grid spec "grid(heg_m1.rmf[cols ENERG_LO" specifies an invalid grid Finished processing grating arfs for ccd_id= 4 5 6 7 ls: acisf03761*HEG_-1.fits: No such file or directory runDmarfadd() ERROR: Could not find grating arfs to concatenate! ...and if you try to run it again, it can't find the parameter file: [cohen@ruby reproc]$ fullgarf Error in parameter file: line 5: field Prompt Could not find the parameter file! (Looked for ) Aborting. ----- Hmmm... looks like maybe two problems - 1. problems finding the fullgarf.par file (though plist works; and you can get it to run by invoking it and then entering params one at a time on the command line) and 2. error when making the grating arfs that maybe has to do with the energy grid - perhaps the more serious, second problem is related to the grating rmfs - try remaking those: http://cxc.harvard.edu/ciao/threads/mkgrmf_aciss/ Remade the heg+1 gRMF - but followed the same directions the same way... output file = heg_p1_rmf.fits Now, trying to make the grating ARF with this rmf (note: "4" is the HEG +1) ---- 11 Aug 05 - with Marc - problem appears to be simply with the setting of the PLIST environment variable... OK, that's fixed (check tg_detect; that should work now too). BUT, there's another problem, with parsing the engrid parameter entry. The issue is the comma separating the ENERG_LO,ENERG_HI - apparently, because it's inside quotes, it is treated as dividing the parameter string in two (i.e. the quotes are, in some sense, being ignored). This isn't the case on Marc's very similar, but not identical, CIAO installation. After trying many things, Marc decided he couldn't fix it. BUT, you can give engrid just the name of the rmf, then when you run fullgarf, it just uses some default column names...and that seems to work. Anyway, for now, I'm just using the garfs Marc made. But this (just giving the rmf file name to engrid and using the default col names) is a klugy workaround. ...while waiting for Marc to test fullgarf now that he thinks he's fixed the PLIST problem, I'll take his garfs (and my rmfs) and do some rudimentary data analysis: Using Victoria's instructions and scripts: Actually, the scripts I have are in the dori/ directory: meg_start.script, analysis.script Q: should I run add_grating_orders? edited two scripts, then ran: [cohen@ruby reproc]$ sherpa meg_start.script ----------------------------------------------------- Welcome to Sherpa: CXC's Modeling and Fitting Program ----------------------------------------------------- Version: CIAO3.2 Type AHELP SHERPA for overview. Type EXIT, QUIT, or BYE to leave the program. Notes: Temporary files for visualization will be written to the directory: /tmp To change this so that these files are not deleted when you exit Sherpa, edit $ASCDS_WORK_PATH in your 'ciao' setup script. Abundances set to Anders & Grevesse The inferred file type is PHA Type II. If this is not what you want, please specify the type explicitly in the data command. Warning: could not find SYS_ERR column WARNING: statistical errors specified in the PHA file. These are currently IGNORED. To use them, type: READ ERRORS "[cols CHANNEL,STAT_ERR]" fitsbin WARNING: backgrounds UP and DOWN are being read from this file, and are being combined into a single background dataset. Warning: could not find SYS_ERR column WARNING: multiple datasets have been input. The next available dataset number is 13. Model parameter prompting is off Analysis Space for Dataset 1: Wavelength Analysis Space for Dataset 2: Wavelength Analysis Space for Dataset 3: Wavelength Analysis Space for Dataset 4: Wavelength Analysis Space for Dataset 5: Wavelength Analysis Space for Dataset 6: Wavelength Analysis Space for Dataset 7: Wavelength Analysis Space for Dataset 8: Wavelength Analysis Space for Dataset 9: Wavelength Analysis Space for Dataset 10: Wavelength Analysis Space for Dataset 11: Wavelength Analysis Space for Dataset 12: Wavelength The inferred file type is ARF. If this is not what you want, please specify the type explicitly in the data command. The inferred file type is ARF. If this is not what you want, please specify the type explicitly in the data command. sherpa> use analysis.script ==> Error bars computed using Chi Gehrels. ==> Error bars computed using Chi Gehrels. did a "fit" then a "lp 2 fit 9 fit 10" then an "uncertainty" Note error message/warning about background. Should we include it? (I say no - source signal is high) *What about doing the analysis in xspec? Looking at Maurice's xspec script - loaddata.xcm copying that in and modifying it Note: no separate background pha2 files... [cohen@ruby reproc]$ xspec Xspec 11.3.1 15:33:32 11-Aug-2005 For documentation, notes, and fixes see http://xspec.gsfc.nasa.gov/ Plot device not set, use "cpd" to set it Type "help" or "?" for further information XSPEC>@loaddata !XSPEC> data 1 acisf03761_000N001_reset_flt1_dstrk_pha2.fits{9}; Net count rate (cts/s) for file 1 8.1031E-02+/- 2.2107E-03 1 data set is in use !XSPEC> data 2 acisf03761_000N001_reset_flt1_dstrk_pha2.fits{10}; Net count rate (cts/s) for file 2 9.1494E-02+/- 2.2459E-03 2 data sets are in use !XSPEC> resp 1 meg_m1.rmf; !XSPEC> tclunknown resp 1 meg_m1.rmf !XSPEC> ::namespace current !XSPEC> response 1 meg_m1.rmf !XSPEC> resp 2 meg_p1.rmf; !XSPEC> tclunknown resp 2 meg_p1.rmf !XSPEC> ::namespace current !XSPEC> response 2 meg_p1.rmf !XSPEC> arf 1 acisf03761MEG_-1_garf.fits; Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) !XSPEC> arf 2 acisf03761MEG_1_garf.fits; Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) XSPEC>cpd /xw & invoke setup.xcm XSPEC>@setup !XSPEC> setplot wave; !XSPEC> query yes; Querying disabled - assuming answer is yes !XSPEC> statistic cstat; Warning : this statistic is only valid for Poisson data XSPEC>ignore 0.0-6.0 ignoring channels 7193 - 8192 in dataset 1 ignoring channels 7193 - 8192 in dataset 2 Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) XSPEC>ignore 6.4-** ignoring channels 1 - 7113 in dataset 1 ignoring channels 1 - 7113 in dataset 2 Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) Note that RESPFILE keyword in ARF is grid(heg_p1_rmf.fits[cols ENERG_LO,ENERG_HI]) XSPEC>plot data XSPEC>model lgauss + pow Model: lgauss<1> + powerlaw<2> Input parameter value, delta, min, bot, top, and max values for ... 18.97 -0.1 1 1 200 200 1:lgauss:waveleng>6.18 0.01 -0.1 0 0 1 1 2:lgauss:sigma>0.1 1 0.01 0 0 1E+24 1E+24 3:lgauss:norm>1.e-3 1 0.01 -3 -2 9 10 4:powerlaw:PhoIndex>2 1 0.01 0 0 1E+24 1E+24 5:powerlaw:norm>1e-2 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: lgauss<1> + powerlaw<2> Model Fit Model Component Parameter Unit Value par par comp 1 1 1 lgauss waveleng A 6.18000 frozen 2 2 1 lgauss sigma A 0.100000 frozen 3 3 1 lgauss norm 1.000000E-03 +/- 0.00000 4 4 2 powerlaw PhoIndex 2.00000 +/- 0.00000 5 5 2 powerlaw norm 1.000000E-02 +/- 0.00000 --------------------------------------------------------------------------- --------------------------------------------------------------------------- C-statistic = 4181.768 using 158 PHA bins. XSPEC>newpar 2 0.01 3 variable fit parameters C-statistic = 5758.117 using 158 PHA bins. Akaike Information Criterion = 6370.030 Bayesian Information Criterion = 6385.584 XSPEC>plot XSPEC>newpar 3 1.e-4 3 variable fit parameters C-statistic = 266.2820 using 158 PHA bins. Akaike Information Criterion = 878.1951 Bayesian Information Criterion = 893.7490 XSPEC>plot XSPEC>thaw 1 Number of variable fit parameters = 4 XSPEC>thaw 2 Number of variable fit parameters = 5 XSPEC>freeze 4 Number of variable fit parameters = 4 XSPEC>fit C-statistic Lvl Fit param # 1 2 3 5 185.657 -2 6.185 3.2047E-03 3.3894E-05 8.9565E-03 185.657 4 6.185 3.2043E-03 3.3894E-05 8.9565E-03 --------------------------------------------------------------------------- Variances and Principal axes : 1 2 3 5 4.89E-11 | 0.00 0.00 -1.00 -0.01 8.72E-08 | 0.00 0.01 -0.01 1.00 1.02E-04 | -1.00 0.00 0.00 0.00 1.73E-02 | 0.00 1.00 0.00 -0.01 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: lgauss<1> + powerlaw<2> Model Fit Model Component Parameter Unit Value par par comp 1 1 1 lgauss waveleng A 6.18543 +/- 0.101031E-01 2 2 1 lgauss sigma A 3.204276E-03 +/- 0.131590 3 3 1 lgauss norm 3.389399E-05 +/- 0.426874E-04 4 4 2 powerlaw PhoIndex 2.00000 frozen 5 5 2 powerlaw norm 8.956515E-03 +/- 0.730482E-03 --------------------------------------------------------------------------- --------------------------------------------------------------------------- C-statistic = 185.6567 using 158 PHA bins. Akaike Information Criterion = 799.5697 Bayesian Information Criterion = 820.3082 XSPEC>plot XSPEC>goodness 100 88.00% of realizations have a fit statistic < 185.7 XSPEC>fit C-statistic Lvl Fit param # 1 2 3 5 182.151 6 6.185 3.2043E-03 3.5703E-05 9.4346E-03 --------------------------------------------------------------------------- Variances and Principal axes : 1 2 3 5 5.42E-11 | 0.00 0.00 -1.00 -0.01 9.21E-08 | 0.00 0.01 -0.01 1.00 1.02E-04 | -1.00 0.00 0.00 0.00 1.73E-02 | 0.00 1.00 0.00 -0.01 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: lgauss<1> + powerlaw<2> Model Fit Model Component Parameter Unit Value par par comp 1 1 1 lgauss waveleng A 6.18543 +/- 0.101029E-01 2 2 1 lgauss sigma A 3.204272E-03 +/- 0.131588 3 3 1 lgauss norm 3.570317E-05 +/- 0.449575E-04 4 4 2 powerlaw PhoIndex 2.00000 frozen 5 5 2 powerlaw norm 9.434590E-03 +/- 0.750541E-03 --------------------------------------------------------------------------- --------------------------------------------------------------------------- C-statistic = 182.1509 using 158 PHA bins. Akaike Information Criterion = 796.0639 Bayesian Information Criterion = 816.8024 XSPEC>goodness 100 86.00% of realizations have a fit statistic < 182.2 After letting the powerlaw index be a free parameter (had to steppar a little), we found our "final" best fit: --------------------------------------------------------------------------- Model: lgauss<1> + powerlaw<2> Model Fit Model Component Parameter Unit Value par par comp 1 1 1 lgauss waveleng A 6.18558 +/- 0.733198E-02 2 2 1 lgauss sigma A 3.688158E-03 +/- 0.928785E-01 3 3 1 lgauss norm 5.050468E-05 +/- 0.450169E-04 4 4 2 powerlaw PhoIndex 2.57994 +/- 1.86366 5 5 2 powerlaw norm 1.339688E-02 +/- 0.169540E-01 --------------------------------------------------------------------------- --------------------------------------------------------------------------- C-statistic = 171.3000 using 158 PHA bins. Akaike Information Criterion = 787.2131 Bayesian Information Criterion = 813.1362 XSPEC>plot XSPEC>iplot PLT> wdata doar21_xspec_618meg.qdp PLT> exit XSPEC>goodness 100 65.00% of realizations have a fit statistic < 171.3 OK, plotted this fit using Kevin's idl code (in ..reproc/idl/) - fit doesn't look great (not narrow fwhm). Before going further, let's make new script loaddata_both.xcm to read in both HEG and MEG data OK... Si XIV fit not so good, when fitting simultaneously. Aside - procedure; discussion with Marc - can use 'uncertainty' in Sherpa, but there's another, steppar like command too...