From a22ed581524e5d39b8d3fed3d0dfe6852d4939be Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Fri, 5 Apr 2024 18:50:27 -0400 Subject: [PATCH 1/9] Add direct loading and thinning of beamstop and TEY signals. --- src/PyHyperScattering/SST1RSoXSDB.py | 70 ++++++++++++++++++++-------- 1 file changed, 50 insertions(+), 20 deletions(-) diff --git a/src/PyHyperScattering/SST1RSoXSDB.py b/src/PyHyperScattering/SST1RSoXSDB.py index 6b1a34cb..e83f1274 100644 --- a/src/PyHyperScattering/SST1RSoXSDB.py +++ b/src/PyHyperScattering/SST1RSoXSDB.py @@ -818,6 +818,7 @@ def loadMonitors( integrate_onto_images: bool = True, useShutterThinning: bool = True, n_thinning_iters: int = 5, + directLoadPulsedMonitors: bool = True ): """Load the monitor streams for entry. @@ -838,51 +839,52 @@ def loadMonitors( useShutterThinning : bool, optional Whether or not to attempt to thin (filter) the raw time streams to remove data collected during shutter opening/closing, by default False As of 9 Feb 2023 at NSLS2 SST1, using useShutterThinning= True for exposure times of < 0.5s is - not recommended because the shutter data is unreliable and too many points will be culled + not recommended because the shutter data is unreliable and too many points will be removed n_thinning_iters : int, optional how many iterations of thinning to perform, by default 5 If the data is becoming too sparse, try fewer iterations + directLoadPulsedMonitors : bool, optional + Whether or not to load the pulsed monitors using direct reading, by default True + This only applies if integrate_onto_images is True; otherwise you'll get very raw data. + If False, the pulsed monitors will be loaded using a shutter-thinning and masking approach as with continuous counters + Returns ------- xr.Dataset xarray dataset containing all monitor streams as data variables mapped against the dimension "time" """ - monitors = None + raw_monitors = None + - # Iterate through the list of streams held by the Bluesky document 'entry' + # Iterate through the list of streams held by the Bluesky document 'entry', and build for stream_name in list(entry.keys()): # Add monitor streams to the output xr.Dataset if "monitor" in stream_name: - if monitors is None: # First one + if raw_monitors is None: # First one # incantation to extract the dataset from the bluesky stream - monitors = entry[stream_name].data.read() + raw_monitors = entry[stream_name].data.read() else: # merge into the to existing output xarray - monitors = xr.merge((monitors, entry[stream_name].data.read())) + raw_monitors = xr.merge((raw_monitors, entry[stream_name].data.read())) # At this stage monitors has dimension time and all streams as data variables # the time dimension inherited all time values from all streams # the data variables (Mesh current, sample current etc.) are all sparse, with lots of nans # if there are no monitors, return an empty xarray Dataset - if monitors is None: + if raw_monitors is None: return xr.Dataset() # For each nan value, replace with the closest value ahead of it in time # For remaining nans, replace with closest value behind it in time - monitors = monitors.ffill("time").bfill("time") + monitors = raw_monitors.ffill("time").bfill("time") # If we need to remap timepoints to match timepoints for data acquisition if integrate_onto_images: try: # Pull out ndarray of 'primary' timepoints (measurement timepoints) - try: - primary_time = entry.primary.data["time"].values - except AttributeError: - if type(entry.primary.data["time"]) == tiled.client.array.DaskArrayClient: - primary_time = entry.primary.data["time"].read().compute() - elif type(entry.primary.data["time"]) == tiled.client.array.ArrayClient: - primary_time = entry.primary.data["time"].read() + primary_time = entry.primary.data["time"].__array__() + primary_time_bins = np.insert(primary_time, 0,0) # If we want to exclude values for when the shutter was opening or closing # This doesn't work for exposure times ~ < 0.5 s, because shutter stream isn't reliable @@ -904,22 +906,50 @@ def loadMonitors( "time" ) + #return monitors # Bin the indexes in 'time' based on the intervales between timepoints in 'primary_time' and evaluate their mean # Then rename the 'time_bin' dimension that results to 'time' monitors = ( - monitors.groupby_bins("time", np.insert(primary_time, 0, 0)) + monitors.groupby_bins("time",primary_time_bins,include_lowest=True) .mean() - .rename_dims({"time_bins": "time"}) + .rename({"time_bins": "time"}) ) - + ''' # Add primary measurement time as a coordinate in monitors that is named 'time' # Remove the coordinate 'time_bins' from the array monitors = ( monitors.assign_coords({"time": primary_time}) .drop_indexes("time_bins") .reset_coords("time_bins", drop=True) - ) - + )''' + + # load direct/pulsed monitors + + for stream_name in list(entry.keys()): + if "monitor" in stream_name and ("Beamstop" in stream_name or "Sample" in stream_name): + # the pulsed monitors we know about are "SAXS Beamstop", "WAXS Beamstop", "Sample Current" + # if others show up here, they could be added + out_name = stream_name.replace("_monitor", "") + mon = entry[stream_name].data.read()[out_name].compute() + SIGNAL_THRESHOLD = 0.1 + threshold = SIGNAL_THRESHOLD*mon.mean('time') + mon_filter = xr.zeros_like(mon) + mon_filter[monthreshold] = 1 + mon_filter.values = scipy.ndimage.binary_erosion(mon_filter) + mon_filtered = mon.where(mon_filter==1) + mon_binned = (mon_filtered.groupby_bins("time",primary_time_bins,include_lowest=True) + .mean() + .rename({"time_bins":"time"}) + ) + + if not directLoadPulsedMonitors: + out_name = 'pl_' + out_name + + monitors[out_name] = mon_binned + monitors = monitors.assign_coords({"time": primary_time}) + + except Exception as e: # raise e # for testing warnings.warn( From 050b4ac9247e05b8cf807d08c8fe5d17988e047a Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Fri, 5 Apr 2024 18:50:41 -0400 Subject: [PATCH 2/9] Update multiindex syntax --- src/PyHyperScattering/SST1RSoXSDB.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/PyHyperScattering/SST1RSoXSDB.py b/src/PyHyperScattering/SST1RSoXSDB.py index e83f1274..1759a71a 100644 --- a/src/PyHyperScattering/SST1RSoXSDB.py +++ b/src/PyHyperScattering/SST1RSoXSDB.py @@ -760,7 +760,7 @@ def subtract_dark(img, pedestal=100, darks=None): monitors = ( monitors.rename({"time": "system"}) .reset_index("system") - .assign_coords(system=index) + .assign_coords(mindex_coords) ) if "system_" in monitors.indexes.keys(): From 4135a3e8788ae00c89f73855649b131a61489a76 Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Fri, 5 Apr 2024 18:51:06 -0400 Subject: [PATCH 3/9] Make baseline checking a bit tolerant and warnings friendlier. --- src/PyHyperScattering/SST1RSoXSDB.py | 31 +++++++--------------------- 1 file changed, 8 insertions(+), 23 deletions(-) diff --git a/src/PyHyperScattering/SST1RSoXSDB.py b/src/PyHyperScattering/SST1RSoXSDB.py index 1759a71a..1275e0d7 100644 --- a/src/PyHyperScattering/SST1RSoXSDB.py +++ b/src/PyHyperScattering/SST1RSoXSDB.py @@ -1066,19 +1066,12 @@ def loadMd(self, run): # print(f'Loading from primary: {phs}, value {primary[rsoxs].values}') except (KeyError, HTTPStatusError): try: - blval = baseline[rsoxs] - if ( - type(blval) == tiled.client.array.ArrayClient - or type(blval) == tiled.client.array.DaskArrayClient - ): - blval = blval.read() + blval = baseline[rsoxs].__array__() md[phs] = blval.mean().round(4) - if blval.var() > 0: + if blval.var() > 1e-4*abs(blval.mean()): warnings.warn( ( - f"While loading {rsoxs} to infill metadata entry for {phs}, found" - f" beginning and end values unequal: {baseline[rsoxs]}. It is" - " possible something is messed up." + f"{phs} changed during scan: {blval}." ), stacklevel=2, ) @@ -1087,28 +1080,20 @@ def loadMd(self, run): md[phs] = primary[md_secondary_lookup[phs]].read() except (KeyError, HTTPStatusError): try: - blval = baseline[md_secondary_lookup[phs]] - if ( - type(blval) == tiled.client.array.ArrayClient - or type(blval) == tiled.client.array.DaskArrayClient - ): - blval = blval.read() + blval = baseline[md_secondary_lookup[phs]].__array__() md[phs] = blval.mean().round(4) - if blval.var() > 0: + if blval.var() > 1e-4*abs(blval.mean()): warnings.warn( ( - f"While loading {md_secondary_lookup[phs]} to infill" - f" metadata entry for {phs}, found beginning and end" - f" values unequal: {baseline[rsoxs]}. It is possible" - " something is messed up." + f"{phs} changed during scan: {blval}." ), stacklevel=2, ) except (KeyError, HTTPStatusError): warnings.warn( ( - f"Could not find {rsoxs} in either baseline or primary. " - f" Needed to infill value {phs}. Setting to None." + f"Could not find {rsoxs} in either baseline or primary while" + f" looking for {phs}. Setting to None." ), stacklevel=2, ) From 985cae1ea7962c30f69238ab1a22a4c25b7ccc99 Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Mon, 8 Apr 2024 04:14:49 -0400 Subject: [PATCH 4/9] Make integrators tolerant of datasets --- .../PFEnergySeriesIntegrator.py | 16 ++++++++++---- src/PyHyperScattering/PFGeneralIntegrator.py | 22 ++++++++++++++----- 2 files changed, 28 insertions(+), 10 deletions(-) diff --git a/src/PyHyperScattering/PFEnergySeriesIntegrator.py b/src/PyHyperScattering/PFEnergySeriesIntegrator.py index f100edf1..99469327 100755 --- a/src/PyHyperScattering/PFEnergySeriesIntegrator.py +++ b/src/PyHyperScattering/PFEnergySeriesIntegrator.py @@ -215,17 +215,25 @@ def integrateImageStack(self,img_stack,method=None,chunksize=None): ''' ''' - + full_stack = None + if isinstance(img_stack,xr.Dataset): + full_stack = img_stack + img_stack = img_stack['I_raw'] if (self.use_chunked_processing and method is None) or method=='dask': func_args = {} if chunksize is not None: func_args['chunksize'] = chunksize - return self.integrateImageStack_dask(img_stack,**func_args) + retval = self.integrateImageStack_dask(img_stack,**func_args) elif (method is None) or method == 'legacy': - return self.integrateImageStack_legacy(img_stack) + retval = self.integrateImageStack_legacy(img_stack) else: raise NotImplementedError(f'unsupported integration method {method}') - + if isinstance(full_stack,xr.Dataset): + retval.name='I_integ' + return xr.merge([full_stack,retval]) + else: + return retval + def createIntegrator(self,en,recreate=False): diff --git a/src/PyHyperScattering/PFGeneralIntegrator.py b/src/PyHyperScattering/PFGeneralIntegrator.py index f7bb0bc1..cecb9fae 100755 --- a/src/PyHyperScattering/PFGeneralIntegrator.py +++ b/src/PyHyperScattering/PFGeneralIntegrator.py @@ -411,18 +411,28 @@ def __init__( def __str__(self): return f"PyFAI general integrator wrapper SDD = {self.dist} m, poni1 = {self.poni1} m, poni2 = {self.poni2} m, rot1 = {self.rot1} rad, rot2 = {self.rot2} rad" - def integrateImageStack(self, img_stack, method=None, chunksize=None): - ''' ''' - - if (self.use_chunked_processing and method is None) or method == 'dask': + def integrateImageStack(self,img_stack,method=None,chunksize=None): + ''' + + ''' + full_stack = None + if isinstance(img_stack,xr.Dataset): + full_stack = img_stack + img_stack = img_stack['I_raw'] + if (self.use_chunked_processing and method is None) or method=='dask': func_args = {} if chunksize is not None: func_args['chunksize'] = chunksize - return self.integrateImageStack_dask(img_stack, **func_args) + retval = self.integrateImageStack_dask(img_stack,**func_args) elif (method is None) or method == 'legacy': - return self.integrateImageStack_legacy(img_stack) + retval = self.integrateImageStack_legacy(img_stack) else: raise NotImplementedError(f'unsupported integration method {method}') + if isinstance(full_stack,xr.Dataset): + retval.name='I_integ' + return xr.merge([full_stack,retval]) + else: + return retval def loadPolyMask(self, maskpoints=[], **kwargs): ''' From 64e41bb3b07ebecd74dee12033829ef80a48a8c2 Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Mon, 8 Apr 2024 04:15:12 -0400 Subject: [PATCH 5/9] First pass at normalizing metadata stream names --- src/PyHyperScattering/SST1RSoXSDB.py | 24 ++++++++++++++++++++++++ 1 file changed, 24 insertions(+) diff --git a/src/PyHyperScattering/SST1RSoXSDB.py b/src/PyHyperScattering/SST1RSoXSDB.py index 1275e0d7..ea2405ad 100644 --- a/src/PyHyperScattering/SST1RSoXSDB.py +++ b/src/PyHyperScattering/SST1RSoXSDB.py @@ -802,6 +802,7 @@ def subtract_dark(img, pedestal=100, darks=None): # retxr = (index,monitors,retxr) monitors.attrs.update(retxr.attrs) retxr = monitors.merge(retxr) + retxr = self._normalize_monitor_names(retxr) if self.use_chunked_loading: # dask and multiindexes are like PEO and PPO. They're kinda the same thing and they don't like each other. @@ -962,6 +963,27 @@ def loadMonitors( ) return monitors + def _normalize_monitor_names(self,run): + """ + Normalize instrument-local monitor names to PyHyper 1.0+ standard names + + """ + rename_dict = { + 'RSoXS Au Mesh Current' : 'i0', + 'NSLS-II Ring Current' : 'src_int' + } + + if run.attrs['rsoxs_config'] == 'saxs': + rename_dict['SAXS Beamstop'] = 'It' + rename_dict['Small Angle CCD Detector_image'] = 'I_raw' + elif run.attrs['rsoxs_config'] == 'waxs': + rename_dict['WAXS Beamstop'] = 'It' + rename_dict['Wide Angle CCD Detector_image'] = 'I_raw' + else: + pass + return run.rename(rename_dict) + + def loadMd(self, run): """ return a dict of metadata entries from the databroker run xarray @@ -1125,6 +1147,7 @@ def loadMd(self, run): md.update(run.metadata) return md + ''' def loadSingleImage(self, filepath, coords=None, return_q=False, **kwargs): """ DO NOT USE @@ -1195,3 +1218,4 @@ def loadSingleImage(self, filepath, coords=None, return_q=False, **kwargs): ) else: return xr.DataArray(img, dims=["pix_x", "pix_y"], attrs=headerdict) + ''' \ No newline at end of file From cb0db2fea416068d5b1cdcc3bb43c6c0e761ca91 Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Mon, 8 Apr 2024 04:16:07 -0400 Subject: [PATCH 6/9] Update minimum xarray version --- requirements.txt | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index d34ce394..df3a01da 100644 --- a/requirements.txt +++ b/requirements.txt @@ -9,7 +9,8 @@ pyfai scikit-image scipy pillow -xarray +# the following minimum version is due to new explicit Coordinate code +xarray>2023.6 tqdm pydata_sphinx_theme # the following pin is due to a security update to numexpr: https://github.com/pydata/numexpr/issues/442 From 03b04299fc3d9b8699562f903c8d0c21bba60db5 Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Tue, 9 Apr 2024 07:42:12 -0400 Subject: [PATCH 7/9] change names some more, standardize on capital I --- src/PyHyperScattering/SST1RSoXSDB.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/PyHyperScattering/SST1RSoXSDB.py b/src/PyHyperScattering/SST1RSoXSDB.py index ea2405ad..36116131 100644 --- a/src/PyHyperScattering/SST1RSoXSDB.py +++ b/src/PyHyperScattering/SST1RSoXSDB.py @@ -969,8 +969,8 @@ def _normalize_monitor_names(self,run): """ rename_dict = { - 'RSoXS Au Mesh Current' : 'i0', - 'NSLS-II Ring Current' : 'src_int' + 'RSoXS Au Mesh Current' : 'I0', + 'NSLS-II Ring Current' : 'Isrc' } if run.attrs['rsoxs_config'] == 'saxs': From 90a00aee9b64ad5726a9bff6bf21ce5319409dfd Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Tue, 9 Apr 2024 07:57:53 -0400 Subject: [PATCH 8/9] Drop Py3.8 support due to xarray dependency --- .github/workflows/main.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index 76434c4d..6f431995 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -12,7 +12,7 @@ jobs: runs-on: ${{ matrix.os}} strategy: matrix: - python-version: [ '3.8', '3.9', '3.10','3.11'] + python-version: [ '3.9', '3.10','3.11', '3.12'] os: [ubuntu-latest, macOS-latest, windows-latest] steps: From b7bc9ce5557e47cc5da35a5140415c20a877f021 Mon Sep 17 00:00:00 2001 From: Peter Beaucage Date: Tue, 9 Apr 2024 08:01:17 -0400 Subject: [PATCH 9/9] Py3.12 not supported yet --- .github/workflows/main.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index 6f431995..a85cf471 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -12,7 +12,7 @@ jobs: runs-on: ${{ matrix.os}} strategy: matrix: - python-version: [ '3.9', '3.10','3.11', '3.12'] + python-version: [ '3.9', '3.10','3.11'] os: [ubuntu-latest, macOS-latest, windows-latest] steps: