example E2 bug fix

docker/mpes/example/E2 Binning of WSe2.ipynb

 %% Cell type:markdown id: tags:
 
 # Binning raw Multidimensional Photoemission Spectroscopy (MPES) data and converting it into the NeXus format
 
 %% Cell type:markdown id: tags:
 
 This example shows how to generate xarray based h5 files from WSe2 trARPES measurement data as detailed in this [paper](https://www.nature.com/articles/s41597-020-00769-8) and how to generate a file in the standardised [MPES NeXus format](https://manual.nexusformat.org/classes/contributed_definitions/NXmpes.html#nxmpes) from it.
 Due to the size of the example file (~6GB) you need at least 40 GB of memory on your computer you're executing this example on. If you just want to have a look on how to convert a pre-binned xarray based h5 file into the NeXus format you may have a look at the simpler [Convert to NeXus example](./E1%20Convert%20to%20NeXus.ipynb), which has lower hardware requirements.
 
 %% Cell type:markdown id: tags:
 
 ## Download RAW data (trARPES data of WSe2)
 
 %% Cell type:markdown id: tags:
 
 Here, we just import `shutil` for extracting the downloaded zip archive and set the main file folder for holding the measurement data.
 
 %% Cell type:code id: tags:
 
 ``` python
 import os
 FDIR = f'{os.getcwd()}/Scan049_1'
 ECAL = f'{os.getcwd()}/energycal_2019_01_08'
 ```
 
 %% Cell type:markdown id: tags:
 
 Since the provided measurement files are rather large (~6GB), they are not directly provided with the example.
 You can [download](https://zenodo.org/record/6369728/files/WSe2.zip) it from zenodo. This may take some time. Place the file in the directory of this notebook afterwards. Under Linux, macOS and in a NORTH container you can directly use the cell below to download the file with curl.
 
 %% Cell type:code id: tags:
 
 ``` python
 ! curl -o WSe2.zip "https://zenodo.org/record/6369728/files/WSe2.zip"
 ```
 
 %% Cell type:markdown id: tags:
 
 Now we extract the measurement files.
 
 %% Cell type:code id: tags:
 
 ``` python
 ! unzip WSe2.zip
 ```
 
 %% Cell type:markdown id: tags:
 
 ## Binning of measurement data
 
 %% Cell type:markdown id: tags:
 
 First we import the necessary packages. For a manual on how to install this dependencies refer to the provided [INSTALL.md](./INSTALL.md) file. If you're running a pre-built docker container or working with the NORTH tools, these dependencies are already available for you.
 
 %% Cell type:code id: tags:
 
 ``` python
 from mpes import base as base, fprocessing as fp, analysis as aly
 import matplotlib.pyplot as plt
 import numpy as np
 import os
 from dask import compute
 import datetime as dt
 import h5py
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Initial data binning for distortion correction
 
 %% Cell type:code id: tags:
 
 ``` python
 parp = fp.parallelHDF5Processor(folder=FDIR)
 parp.gather(identifier=r'/*.h5', file_sorting=True)
 len(parp.files)
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Bin a small range of of files to create a momentum map for distortion correction
 parp.files = parp.files[0:50]
 axes = ['X', 'Y', 't']
 # Important to keep the whole detector area for the initial binning!
 bins = [512, 512, 300]
 ranges = [(0, 2048), (0, 2048), (64000, 68000)]
 parp.parallelBinning(axes=axes, nbins=bins, ranges=ranges, scheduler='threads', ret=False)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Determine correction landmarks
 
 %% Cell type:code id: tags:
 
 ``` python
 # Select an energy slice at the valence band maximum, containing the 6 K-points and the Gamma point as distinct features
 mc = aly.MomentumCorrector(parp.combinedresult['binned'])
 mc.selectSlice2D(slice(165, 175), 2)
 # Extract these high-symmetry points
 mc.featureExtract(mc.slice, sigma=5, fwhm=10, sigma_radius=3)
 mc.view(points=mc.features, annotated=True)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Calculate thin plate spline symmetry correction
 
 %% Cell type:code id: tags:
 
 ``` python
 # Calculate a non-linear coordinate transformation based on thin plate splines that restores 6-fold symmetry
 mc.splineWarpEstimate(image=mc.slice, landmarks=mc.pouter_ord, include_center=True,
                       iterative=False, interp_order=2, update=True)
 mc.view(image=mc.slice_transformed, annotated=True, points={'feats':mc.ptargs}, backend='bokeh', crosshair=True, radii=[75,110,150], crosshair_thickness=0.2)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Image registration
 
 %% Cell type:code id: tags:
 
 ``` python
 # Apply a coordinate translation to move the image into the center of the detector
 mc.coordinateTransform(type='translation', xtrans=70., ytrans=70., keep=True)
 plt.imshow(mc.slice_transformed, origin='lower', cmap='terrain_r')
 plt.axvline(x=256)
 plt.axhline(y=256)
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Rotate the image into a high-symmetry direction
 mc.coordinateTransform( type='rotation', angle=-5, center=(256., 256.), keep=True)
 plt.imshow(mc.slice_transformed, origin='lower', cmap='terrain_r')
 plt.axvline(x=256)
 plt.axhline(y=256)
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Display the final deformation field
 subs = 20
 plt.scatter(mc.cdeform_field[::subs,::subs].ravel(), mc.rdeform_field[::subs,::subs].ravel(), c='b')
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Momentum calibration
 
 %% Cell type:code id: tags:
 
 ``` python
 # Pick one high-symmetry point
 point_b = [252.,255.]
 # Pick the BZ center
 point_a = [308.,346.]
 # give the distance of the two in inverse Angstrom
 distance = np.pi*4/3/3.297
 # Momentum calibration assuming equal scaling along both x and y directions (equiscale=True)
 # Requirements : pixel coordinates of and the momentum space distance between two symmetry points,
 # plus the momentum coordinates
 # of one of the two points
 ext = mc.calibrate(mc.slice_transformed,
                    point_from=point_a,
                    point_to=point_b,
                    dist=distance,
                    equiscale=True,
                    ret=['extent'])
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Display corrected image in calibrated coordinates
 mc.view(image=mc.slice_transformed, imkwds=ext)
 plt.xlabel('$k_x$', fontsize=15)
 plt.ylabel('$k_y$', fontsize=15)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Energy calibration
 
 %% Cell type:code id: tags:
 
 ``` python
 # Bin traces for energy calibration
 axes = ['t']
 bins = [1000]
 ranges = [(63000, 80000)]
 traces, tof = fp.extractEDC(folder=ECAL,
                             axes=axes, bins=bins, ranges=ranges)
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Applied bias voltages (negated, in order to achieve negative binding energies, E-E_F)
 voltages = np.arange(-12.2, -23.2, -1)
 ec = aly.EnergyCalibrator(biases=voltages, traces=traces, tof=tof)
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Normalize traces to maximum
 ec.normalize(smooth=True, span=7, order=1)
 ec.view(traces=ec.traces_normed, xaxis=ec.tof, backend='bokeh')
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Define a TOF feature range, and translate it for each of the traces according to their shift in bias voltage
 rg = [(65000, 65200)]
 ec.addFeatures(traces=ec.traces_normed, refid=0, ranges=rg[0], infer_others=True, mode='append')
 ec.featranges
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # Extract the first peak from each feature range
 ec.featureExtract(traces=ec.traces_normed, ranges=ec.featranges)
 ec.view(traces=ec.traces_normed, peaks=ec.peaks, backend='bokeh')
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Calculate energy calibration
 
 %% Cell type:code id: tags:
 
 ``` python
 # calculate the energy calibration (3rd order polynom). Eref corresponds to the binding energy (E-E_F) of the selected feature in the refid trace.
 refid=5
 Eref=-1.3
 axs = ec.calibrate(ret='all', Eref=Eref, t=ec.tof, refid=refid)
 ec.view(traces=ec.traces_normed, xaxis=ec.calibration['axis'], backend='bokeh')
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Quality of calibration
 
 %% Cell type:code id: tags:
 
 ``` python
 # inspect the quality of the energy calibration
 for i in range(0,len(voltages)):
     plt.plot(ec.calibration['axis']-(voltages[i]-voltages[refid]), ec.traces_normed[i])
 plt.xlim([-15,5])
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Inspect calibration function
 
 %% Cell type:code id: tags:
 
 ``` python
 # energy calibration function vs. TOF
 ec.view(traces=ec.calibration['axis'][None,:], xaxis=ec.tof, backend='matplotlib', show_legend=False)
 plt.scatter(ec.peaks[:,0], ec.biases-ec.biases[refid]+Eref, s=50, c='k')
 plt.xlabel('Time-of-flight', fontsize=15)
 plt.ylabel('Energy (eV)', fontsize=15)
 plt.ylim([-8,6])
 plt.xlim([63400,69800])
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Dataframe processor
 
 %% Cell type:code id: tags:
 
 ``` python
 # create the dask data frame processor
 dfp = fp.dataframeProcessor(datafolder=FDIR)
 dfp.read(source='folder', ftype='h5', timeStamps=True)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Apply energy calibration
 
 %% Cell type:code id: tags:
 
 ``` python
 # apply the energy calibration
 dfp.appendEAxis(E0=ec.calibration['E0'], a=ec.calibration['coeffs'])
 dfp.edf.head(8)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Apply distortion correction
 
 %% Cell type:code id: tags:
 
 ``` python
 # apply the distortion correction
 dfp.applyKCorrection(type='tps_matrix',
                      rdeform_field = mc.rdeform_field,
                      cdeform_field = mc.cdeform_field,
                      X='X', Y='Y', newX='Xm', newY='Ym')
 dfp.edf.head(8)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Apply momentum calibration
 
 %% Cell type:code id: tags:
 
 ``` python
 # apply the momentum calibration
 dfp.appendKAxis(point_b[0], point_b[1], X='Xm', Y='Ym', rstart=parp.binranges[0][0],
                 cstart=parp.binranges[1][0],
                 rstep=parp.binsteps[0],
                 cstep=parp.binsteps[1],
                 fc=mc.calibration['coeffs'][0],
                 fr=mc.calibration['coeffs'][1])
 dfp.edf.head(8)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Apply pump-probe delay axis conversion
 
 %% Cell type:code id: tags:
 
 ``` python
 # calculate the pump-probe delay from the ADC coordinates
 ADCRange = (650, 6900)
 timeRange = (-100, 200)
 dfp.edf['delay'] = timeRange[0] + (dfp.edf['ADC']-ADCRange[0]) *\
     (timeRange[1] - timeRange[0])/(ADCRange[1]-ADCRange[0])
 dfp.edf.head(8)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Bin 4D data in transformed grid
 
 %% Cell type:code id: tags:
 
 ``` python
 # process the 4-dimensional binning
 axes = ['kx', 'ky', 'E', 'delay']
 bins = [50, 50, 100, 21]
 ranges = [(-2, 2), (-2, 2), (-3, 2), (-110, 190)]
 # jittering of energy and ADC should best be done on the bin size of the hardware,
 # not the rebinned bin size. This requires reverse-calculating the jitter amplitudes
 # from the bin sizes.
 TOFrange=[64500,67000]
 e_t_conversion = (base.tof2evpoly(ec.calibration['coeffs'],
                                   ec.calibration['E0'],
                                   TOFrange[0])
                   - base.tof2evpoly(ec.calibration['coeffs'],
                                     ec.calibration['E0'], TOFrange[1])
                  ) / (TOFrange[1] - TOFrange[0])
 d_adc_conversion = (timeRange[1] - timeRange[0]) / (ADCRange[1] - ADCRange[0])
 jitter_amplitude = [0.5,
                     0.5,
                     1*bins[2]/abs(ranges[2][1]-ranges[2][0])*e_t_conversion,
                     1*bins[3]/abs(ranges[3][1]-ranges[3][0])*d_adc_conversion]
 dfp.distributedBinning(axes=axes,
                        nbins=bins,
                        ranges=ranges,
                        scheduler='threads',
                        ret=False,
                        jittered=True,
                        jitter_amplitude=jitter_amplitude)
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Create metatada structure
 This adds additional metadata to the h5 file. This data may also be provided through additional ELN entries through a NOMAD instance or with a handwriten file directly to the mpes parser.
 
 %% Cell type:code id: tags:
 
 ``` python
 metadata = {}
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 # manual Meta data
 #General
 metadata['experiment_summary'] = 'WSe2 XUV NIR pump probe data.'
 metadata['entry_title'] = 'Valence Band Dynamics - 800 nm linear s-polarized pump, 0.6 mJ/cm2 absorbed fluence'
 metadata['experiment_title'] = 'Valence band dynamics of 2H-WSe2'
 
 #User
 # Fill general parameters of NXuser
 # TODO: discuss how to deal with multiple users?
 metadata['user0'] = {}
 metadata['user0']['name'] = 'Julian Maklar'
 metadata['user0']['role'] = 'Principal Investigator'
 metadata['user0']['affiliation'] = 'Fritz Haber Institute of the Max Planck Society'
 metadata['user0']['address'] = 'Faradayweg 4-6, 14195 Berlin'
 metadata['user0']['email'] = 'maklar@fhi-berlin.mpg.de'
 
 #NXinstrument
 metadata['instrument'] = {}
 #analyzer
 metadata['instrument']['analyzer']={}
 metadata['instrument']['analyzer']['slow_axes'] = "delay" # the scanned axes
 metadata['instrument']['analyzer']['spatial_resolution'] = 10.
 metadata['instrument']['analyzer']['energy_resolution'] = 110.
 metadata['instrument']['analyzer']['momentum_resolution'] = 0.08
 metadata['instrument']['analyzer']['working_distance'] = 4.
 metadata['instrument']['analyzer']['lens_mode'] = "6kV_kmodem4.0_30VTOF.sav"
 
 #probe beam
 metadata['instrument']['beam']={}
 metadata['instrument']['beam']['probe']={}
 metadata['instrument']['beam']['probe']['incident_energy'] = 21.7
 metadata['instrument']['beam']['probe']['incident_energy_spread'] = 0.11
 metadata['instrument']['beam']['probe']['pulse_duration'] = 20.
 metadata['instrument']['beam']['probe']['frequency'] = 500.
 metadata['instrument']['beam']['probe']['incident_polarization'] = [1, 1, 0, 0] # p pol Stokes vector
 metadata['instrument']['beam']['probe']['extent'] = [80., 80.]
 #pump beam
 metadata['instrument']['beam']['pump']={}
 metadata['instrument']['beam']['pump']['incident_energy'] = 1.55
 metadata['instrument']['beam']['pump']['incident_energy_spread'] = 0.08
 metadata['instrument']['beam']['pump']['pulse_duration'] = 35.
 metadata['instrument']['beam']['pump']['frequency'] = 500.
 metadata['instrument']['beam']['pump']['incident_polarization'] = [1, -1, 0, 0] # s pol Stokes vector
 metadata['instrument']['beam']['pump']['incident_wavelength'] = 800.
 metadata['instrument']['beam']['pump']['average_power'] = 300.
 metadata['instrument']['beam']['pump']['pulse_energy'] = metadata['instrument']['beam']['pump']['average_power']/metadata['instrument']['beam']['pump']['frequency']#µJ
 metadata['instrument']['beam']['pump']['extent'] = [230., 265.]
 metadata['instrument']['beam']['pump']['fluence'] = 0.15
 
 #sample
 metadata['sample']={}
 metadata['sample']['preparation_date'] = '2019-01-13T10:00:00+00:00'
 metadata['sample']['sample_history'] = 'Cleaved'
 metadata['sample']['chemical_formula'] = 'WSe2'
 metadata['sample']['description'] = 'Sample'
 metadata['sample']['name'] = 'WSe2 Single Crystal'
 
 metadata['file'] = {}
 metadata['file']["trARPES:Carving:TEMP_RBV"] = 300.
 metadata['file']["trARPES:XGS600:PressureAC:P_RD"] = 5.e-11
 
 #Sample positions
 metadata['file']["trARPES:Carving:TRX.RBV"] = np.nan
 metadata['file']["trARPES:Carving:TRY.RBV"] = np.nan
 metadata['file']["trARPES:Carving:TRZ.RBV"] = np.nan
 metadata['file']["trARPES:Carving:THT.RBV"] = np.nan
 metadata['file']["trARPES:Carving:PHI.RBV"] = np.nan
 metadata['file']["trARPES:Carving:OMG.RBV"] = np.nan
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Run the following cell to store metadata from EPICS archive only if outside the FHI network
 
 %% Cell type:code id: tags:
 
 ``` python
-metadata['file'] = {}
 metadata['file']["KTOF:Lens:Extr:I"] = -0.12877
 metadata['file']["KTOF:Lens:UDLD:V"] = 399.99905
 metadata['file']["KTOF:Lens:Sample:V"] = 17.19976
 metadata['file']["KTOF:Apertures:m1.RBV"] = 3.729931
 metadata['file']["KTOF:Apertures:m2.RBV"] = -5.200078
 metadata['file']["KTOF:Apertures:m3.RBV"] = -11.000425
 ```
 
 %% Cell type:markdown id: tags:
 
 ### Generate xarray
 
 %% Cell type:code id: tags:
 
 ``` python
 res_xarray = dfp.gather_metadata(metadata_dict=metadata.copy(), ec=ec, mc=mc)
 ```
 
 %% Cell type:markdown id: tags:
 
 ## Create a NeXus file from a xarray
 This conversion basically follows the same procedure as in the [convert to NeXus example](./E1%20Convert%20to%20Nexus.ipynb). Please refer to this notebook for details on the convert function. Here, we are using the objects keywords of `convert` to pass the generated xarray directly, instead of loading a h5 datafile.
 
 %% Cell type:code id: tags:
 
 ``` python
 from nexusparser.tools.dataconverter.convert import convert
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 convert(input_file=["config_file.json"],
         objects=res_xarray,
         reader='mpes',
         nxdl='NXmpes',
         output='WSe2.mpes.nxs')
 ```
 
 %% Cell type:markdown id: tags:
 
 ## View the data with H5Web
 H5Web is a tool for visualizing any data in the h5 data format. Since the NeXus format builds opon h5 it can be used to view this data as well. We just import the package and call H5Web with the output filename from the convert command above. For an analysis on NeXus data files please refer to [analysis example](./E3%20pyARPES%20analysis.ipynb).
 
 You can also view this data with the H5Viewer or other tools from your local filesystem.
 
 %% Cell type:code id: tags:
 
 ``` python
 from jupyterlab_h5web import H5Web
 ```
 
 %% Cell type:code id: tags:
 
 ``` python
 H5Web('WSe2.mpes.nxs')
 ```