################################################## ## Herschel Data Processing Workshop - Feb 2012 ## ## ## ## HIFI Hands-on sessions - 23/24 Feb. 2012 ## ## ## ## DT, CmC, 22-02-2012 ## ################################################## ################################################## ## 1. DATA ACCESS ## ################################################## ##1.1 Single point Obsid: 1342190183, DBS Fast, CO 5-4 in LDN1577 #1.1.1 Access data: option 1 - get from HSA obspoint_hipe6 = getObservation(1342190183,useHsa=True) #Data as processed with the latest bulk reprocessing #1.1.2. Access data: option 2 - get from pool. ##Your data need to be in ~/.hcss/lstore obspoint_hipe6 = getObservation(obsid = 1342190183, poolName='1342190183_hipe6') #Data are here processed with HIPE 6.1.0 ##Inspect data - USE OF GUI ## - check science data ## - quick-look: overplot H and V ## - quick-look: overplot WBS and HRS ## - manual flag inspection ##1.2 Spectral scan Obsid: 1342191505, DBS, OrionS mini-scan in 1a #1.2.1 Access data: option 1 - get from HSA obsscan_hipe6 = getObservation(1342191505,useHsa=True) #Data as processed with the latest bulk reprocessing #1.2.2 Access data: option 2 - get from pool. ##Your data need to be in ~/.hcss/lstore obsscan_hipe6 = getObservation(obsid = 1342191505, poolName='1342191505_hipe6') #Data are here processed with HIPE 6.1.0 ##Inspect data - use GUI ## - check science data ##1.3 Mapping Obsid: 1342210099, OTF map in IRC+10216, CO 5-4 ## 1342205481, DBS raster in C+ on S140 #1.3.1 Access data: option 1 - get from HSA obsmap_hipe6 = getObservation(1342210099,useHsa=True) #Data as processed with the latest bulk reprocessing #1.3.2 Access data: option 2 - get from pool. ##Your data need to be in ~/.hcss/lstore obsmap_hipe6 = getObservation(obsid = 1342210099, poolName='1342210099_hipe6') #Data are here processed with HIPE 6.1.0 ##Inspect data - use GUI ## - check science data ## - inspect spectra in a mosaic display ## - inspect cubes ################################################## ## 2. DATA REPROCESSING ## ################################################## #################################### ##2.1 Single point Obsid: 1342190183 ##Click on the variable obspoint_hipe6 and open ##the hifiPipeline GUI in the Applicable tasks ##Reprocess data: reprocess only WBS data #2.1.1 Default parameters, no calibration update #The equivalent command-line call is: obspoint_reprocess = hifiPipeline(obs=obspoint_hipe6, apids=["WBS-H"], fromLevel=0, \ upToLevel=2.0) #Inspect data: # - note the creation of browse product #Save the newly processed data, here together with calibration files #Option 1: no pool specified: the data will be saved in ~/.hcss/lstore/'obsid' saveObservation(obspoint_reprocess, saveCalTree=True) #Option 2: you can specify a new or an existing pool saveObservation(obspoint_reprocess, poolName = '1342190183_reprocessed_hipe8', saveCalTree=True) ##Read back the observation: same as point 2 above. #2.1.2 Use a different calibration #First define where to fetch the latest calibration files cal = configureHifiPipeline(useHsa=True) #The equivalent command-line call is: obspoint_reprocess = hifiPipeline(obs=obspoint_hipe6, apids=["WBS-H"], fromLevel=0, \ upToLevel=2.0,cal=True, palStore=cal, calVersion="HIFI_CAL_7_0") #Inspect data: # - note the addition of calVersion #2.1.3. Default parameter, but configure pipeline algorithm ##For example: run up to level2 but do not average level2 products ##Easier from the GUI, but the equivalent command would be: #A. obtain a params object htp = obspoint_hipe6.level.get("level0").getProduct("WBS-H") params = herschel.hifi.pipeline.generic.PipelineConfiguration.getConfig(htp) #B. configure it params.setParameter("doAvg","ignore",True) #C. run the pipeline while passing the params object to it: obspoint_reprocess = hifiPipeline(obs=obspoint_hipe6, apids=["WBS-H"], fromLevel=0, \ upToLevel=2.0, params=params) ###################################### ##2.2 Spectral scan Obsid: 1342191505 #2.2.1 Reprocess with alternative calibration scheme ##Click on the variable obsscan_hipe6 and open ##the hifiPipeline GUI in the Applicable tasks # First make a copy of an unreprocessed dataset in order to compare, here we # take the first scan in the level 2 WBS-H-USB ds1 = obsscan_hipe6.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product["0002"].copy() # # obtain a params object htp = obsscan_hipe6.level.get("level0").getProduct("WBS-H") params = herschel.hifi.pipeline.generic.PipelineConfiguration.getConfig(htp) # configure it params.setParameter("doFilterLoads","ignore",False) params.setParameter("doFilterLoads","filterMethod","cubic_splines") # run the pipeline while passing the params object to it: obsscan_reprocess = hifiPipeline(obs=obsscan_hipe6, apids=["WBS-H"], fromLevel=0, \ upToLevel=2.0, params=params) ds2 = obsscan_reprocess.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product["0002"].copy() # #Compare new data with old data by plotting the first dataset in the level 2 WBS-H-USB # and dragging the ds1 created above into the same plot - check that the standing wave has disminished ############################### ##2.3 Mapping Obsid: 1342210099 #Repipeline to get new cubes from herschel.hifi.pipeline.level0 import DoUplink doUplink=DoUplink() backend = "WBS-H" htp=obsmap_hipe6.refs["level1"].product.refs[backend].product htp1=doUplink(htp=htp.copy(), auxiliary=obsmap_hipe6.getAuxiliary()) l1=obsmap_hipe6.refs["level1"].product l1.setProduct(backend,htp1) obsmap_hipe6.setProduct("level1",l1) obsmap_reprocess = hifiPipeline(obs=obsmap_hipe6, fromLevel=1, upToLevel=2,apids=["WBS-H"]) ################################################## ## 3. DATA MANIPULATION ## ################################################## ################################################## ## 3.1 FRINGE REMOVAL AND FLAG HANDLING ## ################################################## ###################################### ##3.1.1 Single point Obsid: 1342190183 ##Example in Single point data, use data reprocessed with HIPE 8 ##Here, the mask for line will be assign by us obspoint_hipe8 = getObservation(obsid = 1342190183, poolName='1342190183_hipe8') ##Open the FitHifiFringe GUI - e.g.: click on 'obspoint_hipe8' in the ##VARIABLE column and look at the Applicable task in the "Tasks" column ##3.1.1.1: With mask window indicated in the GUI but no baseline removal. Concentrate on WBS-H-USB ##Here we indicate e.g. a line window between 576.2 and 576.4 GHz ##The equivalent command line will be: obspoint_defringed = fitHifiFringe(obs_or_htp=obspoint_hipe8, product='WBS-H-USB', nfringes=1, usermask=(576.2, 576.4), automask=False) ##Repeat for WBS-V and save baseline-subtracted data obspoint_defringed = fitHifiFringe(obs_or_htp=obspoint_hipe8, product='WBS-V-USB', nfringes=1, usermask=(576.2, 576.4), automask=False) saveObservation(obspoint_defringed, poolName='1342190183_defringed',saveCalTree=True) ##3.1.1.2: With mask window indicated by BRIGHT_LINE flag but no baseline removal. Concentrate on WBS-H-USB ##First assign pixel flags: ## - you can do it in the GUI - need to know the flag value: 536870912 = 2^29 for BRIGHT_LINE (see documentation) ## - otherwise the command line needs to know the sub-band nb and channel indices for the window ds_point = obspoint_hipe8.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product["0001"] flagPixels(ds=ds_point, selection=[0], mask={2:range(681,1081)}, flag = 536870912) ##Open ds_point in the Spectrum explorer and display the flags (click on the "Flag" icon). ##The selected window should show up in colour. ##Then run defringing and observe how the flagged range is masked obspoint_defringed = fitHifiFringe(sds1=ds_point, nfringes=1, automask=False) ####################################### ##3.1.2 Spectral Scan Obsid: 1342191505 ##Example in Spectral scan data reprocessed with HIPE 8 obsscan_hipe8 = getObservation(obsid = 1342191505, poolName='1342191505_hipe8') ##Open the FitHifiFringe GUI - e.g.: click on 'obsscan_hipe8' in the ##VARIABLE column and look at the Applicable task in the "Tasks" column ##This time, there are plenty of lines, so we will let the task choose ##mask at each LO setting ##3.1.2.1: with automask but no baseline removal. Concentrate on WBS-H-USB ##The equivalent command line will be: obsscan_defringed = fitHifiFringe(obs_or_htp=obsscan_hipe8, nfringes=2) saveObservation(obsscan_defringed, poolName='1342191505_defringed',saveCalTree=True) ##3.1.2.2: with automask and baseline removal. Concentrate on WBS-H-USB ##The equivalent command line will be: obsscan_defringed = fitHifiFringe(obs_or_htp=obsscan_hipe8, nfringes=2, sub_base=True) ##Compare the output of obsscan_defringed with that of obsscan_reprocess from step 2.2.1 above ################################################## ## 3.2 BASELINE CORRECTION ## ################################################## ####################################### ##3.2.1 Spectral Scan Obsid: 1342191505 ##Example in Spectral scan data, after fringe removal obsscan_defringed = getObservation('1342191505', poolName='1342191505_defringed') ##We are concentrating on the WBS-H data htp_wbsh_scan = obsscan_defringed.refs["level2"].product.refs["WBS-H-USB"].product.copy() ##In this spectral scan, there are plenty of lines so it is handy to have the task identify the masks. ##The masks will propagate from one LO setting to another ##Note how the LSB will be done automatically, based on the choice made for the USB fitBaseline(obs=obsscan_defringed, htp=htp_wbsh_scan, backends=["WBS-H"], order=3, domask = 2, interactive=True) htp_baselined = fitBaseline.result_htp obsscan_baselined = fitBaseline.result_obs ##Save baselined data in order to deconvolve at later stage saveObservation(obsscan_baselined, poolName='1342191505_baselined',saveCalTree=True) ######################################################################## ##Additional example for illustration only: case of broad lines ##M82 single point: broad line, whereby the automask does not work well ##obs = getObservation(1342207319,useHsa=True) ######################################################################## ################################################## ## 3.3 DECONVOLUTION ## ################################################## ####################################### ##3.3.1 Spectral Scan Obsid: 1342191505 ##We will work from the defringed, baselined, data obsscan_baselined = getObservation('1342191505', poolName='1342191505_baselined') ##Open the doDeconvolution GUI - e.g.: click on 'obsscan_baselined' in the ##VARIABLE column and look at the Applicable task in the "Tasks" column #The equivalent command line is: obsscan_deconvolved = doDeconvolution(obs=obsscan_baselined, polarization="H_POL", plot_dsb = 'PLOT_ALL') ######################################################################## ##Additional example for illustration only: case with lots of spurs ##obs = getObservation(1342231503,useHsa=True) ######################################################################## ################################################## ## 3.4 SPECTRAL TOOLS: SOME EXAMPLES ## ################################################## ####################################### ##3.4.1 Spectral Scan Obsid: 1342191505 #Open the deconvolution output in the Spectrum explorer decon_spec = obsscan_deconvolved["ssb"].copy() ##3.4.1.1 Smooth the data #The equivalent command line is: decon_spec_smoothed = smooth(ds=decon_spec, filter='Gaussian', width=4.0, edge='ZEROES', variant='flag-weight', overwrite=False) ##3.4.1.2 Convert into km/s - reference velocity taken at rest frequency of CS 10-9 at 489.750921 GHz #The equivalent command line is: decon_spec_velo = convertWavescale(ds=decon_spec_smoothed, to='km s-1', reference=489.750921, referenceUnit='GHz', overwrite=False) ##3.4.1.3 Convert into cm-1 #The equivalent command line is: decon_spec_wavenumber = convertWavescale(ds=decon_spec_smoothed, to='cm-1', referenceUnit='GHz', overwrite=False) ###################################### ##3.4.2 Simple point Obsid: 1342190183 ##3.4.2.1 Convert to Jy: example on single spectrum data obspoint_defringed = getObservation(obsid = 1342190183, poolName='1342190183_defringed') ds1 = obspoint_defringed.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product["0001"].copy() ##Some telescope efficiencies are needed - the best is to ##pass them from the calibration product ##You also need to make an assumption about your source size cal=obspoint_defringed.calibration ds1_Jy = convertK2Jy(ds=ds1,cal=cal, size=0) ##3.4.2.2 Average the two polarisations ds1 = obspoint_defringed.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product["0001"].copy() ds2 = obspoint_defringed.refs["level2"].product.refs["WBS-V-USB"].product.refs["box_001"].product["0001"].copy() ds_average = polarPair(ds1=ds1, ds2=ds2) ################################################## ## 3.5 CUBE MANIPULATION ## ################################################## ###################################### ##3.5.1 Spectral Map Obsid: 1342210099 ##3.5.1.1 Regridding of cleaned data ##We will use baseline-corrected data from IRC+10216 provided in one of the pools obsmap_baselined = getObservation('1342210099', poolName='1342210099_baselined') ##Gridding will be done per Hifi Timeline product htp_wbsh_map = obsmap_baselined.refs["level2"].product.refs["WBS-H-USB"].product ##Open doGridding GUI: optimum parameters are open by default ##We will limit ourselves to sub-band 1 where the strong CO 5-4 is located. #The equivalent command line is: doGridding(htp=htp_wbsh_map, subbands=[1], beam=[37.12009380988368,37.12009380988368], \ pixelSize=[15.4667057541182,15.4667057541182], weightMode="equal", filterType="gaussian", \ xFilterParams=[2.1600000027298036,0.43240406770959333 ,2.0], yFilterParams=[2.1600000027298036,0.43240406770959333,2.0], \ datasetTypes=["science"], channels=[[0,2266]], refPixelCoordinates=[146.9886228165384,13.278050038400679]) cube_regridded = doGridding.cube ##3.5.1.2: Open cube with the Cube Analysis Tool Box GUI ##Some examples: ## - let's create an integrated map intensity of the CO (5-4) line ## - let's create a collapsed spectrum over the whole cube ## - let's create a pv-diagram with reference frequency the CO 5-4 rest freq. at 576.26793 GHz ################################################## ## 3.6 DATA FITTING ## ################################################## ####################################### ##3.6.1 Spectral scan Obsid: 1342191505 obsscan_baselined = getObservation('1342191505', poolName='1342191505_baselined') ds_sscan_8 = obsscan_baselined.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product["0008"] ##Open spectrum in the Spectral Explorer and open the Spectral toolbox ##Go to SpectrumFitterGUI - the rest will be done in a GUI fashion ################################################## ## 4 DATA EXPORT ## ################################################## #################################### ##4.1 Single Point Obsid: 1342190183 ##4.1.1 Export to GILDAS/Class ##Here for example, export both H and V USB products obspoint_defringed = getObservation(obsid = 1342190183, poolName='1342190183_defringed') htp_husb = obspoint_defringed.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product htp_vusb = obspoint_defringed.refs["level2"].product.refs["WBS-H-USB"].product.refs["box_001"].product hclass = hiClass() hclass.add(htp_husb) hclass.add(htp_vusb) ##Change export directory according to your needs hclass.saveToFits("/Users/dteyssier/1342190183_hipe8_USB_Class.fits") ##Now you need to import the data into Class - see documentation (chapter 18.3) ##4.1.2 Export to standard FITS ##Change export directory according to your needs simpleFitsWriter(product=htp_husb,file="/Users/dteyssier/1342190183_hipe8_WBS_H_USB.fits") simpleFitsWriter(product=htp_vusb,file="/Users/dteyssier/1342190183_hipe8_WBS_V_USB.fits")