18jun12 in xspec/2012/zpup_rgs/grouping/ (but older sas11 version, directory may move) the upshot is it looks like Vero's code is tuned to the structure of ciao pha files (where, e.g., there are columns for the low and high wavelengths of each bin; whereas for xmm pha files, the wavelength solution appears to be encoded analytically in the header). Aside from this, it looks like the first thing to do is to run mrdfits to read the pha file into a structure: IDL> data = mrdfits("acisf00638_000N003_repro_pha2.fits",1,head) IDL> help,data,/struct ** Structure <87a9d24>, 13 tags, length=229392, data length=229392, refs=1: SPEC_NUM INT 1 TG_M INT -3 TG_PART INT 1 TG_SRCID INT 1 X FLOAT 4070.49 Y FLOAT 4162.05 CHANNEL INT Array[8192] COUNTS INT Array[8192] STAT_ERR FLOAT Array[8192] BACKGROUND_UP INT Array[8192] BACKGROUND_DOWN INT Array[8192] BIN_LO DOUBLE Array[8192] BIN_HI DOUBLE Array[8192] IDL> data = mrdfits("r1m1_grp_min40.pha",1,head) MRDFITS: Binary table. 6 columns by 3600 rows. IDL> help,data,/struct ** Structure <87a7a7c>, 6 tags, length=20, data length=20, refs=1: CHANNEL LONG 1 COUNTS LONG 0 QUALITY INT 1 GROUPING INT 1 AREASCAL FLOAT 0.00000 BACKSCAL FLOAT 0.00000 Here the help.../struct functionality shows that chandra and xmm pha files have pretty different column structures Note also that individal variables in a structure can be accessed using structure.variable syntax: IDL> print, data.grouping 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 -1 -1 -1 -1 -1 -1 -1 -1 -1 1 -1 -1 -1 -1 -1 1 1 -1 1 1 -1 1 1 -1 1 -1 1 -1 1 1 -1 1 1 1 -1 1 -1 1 -1 1 -1 1 1 1 -1 1 1 -1 1 1 1 1 1 -1 1 -1 1 -1 1 -1 1 -1 1 1 -1 1 1 1 -1 1 -1 1 1 1 1 1 1 1 -1 1 1 1 1 1 -1 1 -1 1 1 -1 1 1 1 1 1 1 1 1 Here we see that for the RGS data file, the grouping column (made with grppha) is just 1 and -1 (showing group boundaries?) ---------- Now, evolving Vero's code... or really just starting from scratch and simply flipping "1" to "-1" for some user-determined bins, and then write out a new fits file. IDL> data = mrdfits("rgs1m1_combined_he_bin.pha",1,head) IDL> help,data,/struct ** Structure <9d40494>, 6 tags, length=20, data length=20, refs=1: CHANNEL LONG 1 COUNTS LONG 0 QUALITY INT 1 GROUPING INT 1 AREASCAL FLOAT 0.00000 BACKSCAL FLOAT 0.00000 IDL> data_new = data IDL> help,data_new,/struct ** Structure <9d40494>, 6 tags, length=20, data length=20, refs=2: CHANNEL LONG 1 COUNTS LONG 0 QUALITY INT 1 GROUPING INT 1 AREASCAL FLOAT 0.00000 BACKSCAL FLOAT 0.00000 SO, I can copy a structure IDL> mwrfits,data_new,'test.fits',head,/create and write it out to a new fits file (I checked the structure and spot checked the grouping) NEXT: change the first and second 1 to a -1 IDL> data = mrdfits("rgs1m1_combined_he_bin.pha",1,head) % Compiled module: MRDFITS. % Compiled module: FXPOSIT. % Compiled module: FXMOVE. % Compiled module: MRD_HREAD. % Compiled module: FXPAR. % Compiled module: GETTOK. % Compiled module: VALID_NUM. % Compiled module: MRD_SKIP. % Compiled module: MATCH. % Compiled module: MRD_STRUCT. MRDFITS: Binary table. 6 columns by 3600 rows. % Compiled module: IS_IEEE_BIG. IDL> data_new = data IDL> data_new[1:2].grouping = -1 IDL> mwrfits,data_new,'test2.fits',head,/create % Compiled module: MWRFITS. % Compiled module: FXADDPAR. % Compiled module: FXPARPOS. % Compiled module: DETABIFY. % Compiled module: HOST_TO_IEEE. IDL> exit Yes, it worked: IDL> data = mrdfits("test2.fits",1,head) % Compiled module: MRDFITS. % Compiled module: FXPOSIT. % Compiled module: FXMOVE. % Compiled module: MRD_HREAD. % Compiled module: FXPAR. % Compiled module: GETTOK. % Compiled module: VALID_NUM. % Compiled module: MRD_SKIP. % Compiled module: MATCH. % Compiled module: MRD_STRUCT. MRDFITS: Binary table. 6 columns by 3600 rows. % Compiled module: IS_IEEE_BIG. IDL> print,data.grouping 1 -1 -1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 --------- 23jun12 -- TESTING IT OUT ON GROUP 40 PHA IDL> data = mrdfits("r1m1_grp_min40.pha",1,head) IDL> data_new = data IDL> data_new[2450:2470].grouping = -1 IDL> mwrfits,data_new,'test3.fits',head,/create That range is what Jake found for zeta Pup RGS1 N VI (with the gap shift applied) NOW UNPACKING THE SAS12 REDUCTION (AND STARTING OVER WITH THE GRPPHA GROUPING AND THEN THE GROUPING BY HAND) [cohen@hven zpup_rgs]$ ftools6.12 [cohen@hven zpup_rgs]$ grppha Please enter PHA filename[r1m1.pha] r1m1.pha Please enter output filename[r1m1_grp_OandN.pha] r1m1_min40.pha GRPPHA[group 1730 1850 121 2420 2570 151] group min 40 GRPPHA[exit] ... written the PHA data Extension ...... exiting, changes written to file : r1m1_min40.pha ** grppha 3.0.1 completed successfully DID THE SAME FOR RGS2 NOW DO THE HE-LIKE BINS BY HAND IDL> data = mrdfits("r1m1_min40.pha",1,head) IDL> data_new = data IDL> data_new[105:118].grouping = -1 IDL> data_new[300:338].grouping = -1 IDL> data_new[495:560].grouping = -1 IDL> data_new[930:995].grouping = -1 IDL> data_new[1740:1840].grouping = -1 IDL> data_new[2450:2620].grouping = -1 IDL> mwrfits,data_new,'r1m1_min40_binHe.pha',head,/create AND FOR RGS2 IDL> data = mrdfits("r2m1_min40.pha",1,head) IDL> data_new = data IDL> data_new[251:287].grouping = -1 IDL> data_new[496:561].grouping = -1 IDL> data_new[791:851].grouping = -1 IDL> data_new[2454:2624].grouping = -1 IDL> mwrfits,data_new,'r2m1_min40_binHe.pha',head,/create !NOTE THAT THE ABOVE RANGES HAVE BEEN UPDATED A FEW TIMES, AND ARE CURRENT AS OF 4JUL12