Skip to content

Commit bb22360

Browse files
aglowackiaglowacki
authored
Merge pull request #195 from aglowacki/master
Refactor snip
2 parents a457ec6 + d71a5c7 commit bb22360

File tree

3 files changed

+116
-40
lines changed

3 files changed

+116
-40
lines changed

‎src/data_struct/spectra.h‎

Lines changed: 53 additions & 39 deletions
Original file line numberDiff line numberDiff line change
@@ -52,6 +52,7 @@ POSSIBILITY OF SUCH DAMAGE.
5252
#include "core/defines.h"
5353
#include <Eigen/Core>
5454
#include <vector>
55+
#include <algorithm>
5556
#include <functional>
5657

5758
namespace data_struct
@@ -254,6 +255,7 @@ class Spectra : public ArrayTr<_T>
254255
TEMPLATE_CLASS_DLL_EXPORT Spectra<float>;
255256
TEMPLATE_CLASS_DLL_EXPORT Spectra<double>;
256257

258+
//-----------------------------------------------------------------------------
257259
// ----------------------------------------------------------------------------
258260

259261
template<typename T_real>
@@ -345,31 +347,17 @@ DLL_EXPORT ArrayTr<T_real> snip_background(const Spectra<T_real> * const spectra
345347

346348
// FIRST SNIPPING
347349
int no_iterations = 2;
348-
349-
int max_of_xmin = (std::max)(xmin, (T_real)0.0);
350-
int min_of_xmax = (std::min)(xmax, T_real(spectra->size() - 1));
350+
351+
int lo_limit = std::max<int>(xmin, 0);
352+
int up_limit = std::min<int>(xmax, (spectra->size() - 1));
351353
for (int j = 0; j < no_iterations; j++)
352354
{
353355
for (long int k = 0; k < background.size(); k++)
354356
{
355-
long int lo_index = k - current_width[k];
356-
long int hi_index = k + current_width[k];
357-
if (lo_index < max_of_xmin)
358-
{
359-
lo_index = max_of_xmin;
360-
}
361-
if (lo_index > min_of_xmax)
362-
{
363-
lo_index = min_of_xmax;
364-
}
365-
if (hi_index > min_of_xmax)
366-
{
367-
hi_index = min_of_xmax;
368-
}
369-
if (hi_index < max_of_xmin)
370-
{
371-
hi_index = max_of_xmin;
372-
}
357+
int lo_index = std::max<int>(k - current_width[k], lo_limit);
358+
lo_index = std::min<int>(lo_index, up_limit);
359+
int hi_index = std::max<int>(k + current_width[k], lo_limit);
360+
hi_index = std::min<int>(k + current_width[k], up_limit);
373361
T_real temp = (background[lo_index] + background[hi_index]) / (T_real)2.0;
374362
if (background[k] > temp)
375363
{
@@ -382,24 +370,10 @@ DLL_EXPORT ArrayTr<T_real> snip_background(const Spectra<T_real> * const spectra
382370
{
383371
for (long int k = 0; k < background.size(); k++)
384372
{
385-
long int lo_index = k - current_width[k];
386-
long int hi_index = k + current_width[k];
387-
if (lo_index < max_of_xmin)
388-
{
389-
lo_index = max_of_xmin;
390-
}
391-
if (lo_index > min_of_xmax)
392-
{
393-
lo_index = min_of_xmax;
394-
}
395-
if (hi_index > min_of_xmax)
396-
{
397-
hi_index = min_of_xmax;
398-
}
399-
if (hi_index < max_of_xmin)
400-
{
401-
hi_index = max_of_xmin;
402-
}
373+
int lo_index = std::max<int>(k - current_width[k], lo_limit);
374+
lo_index = std::min<int>(lo_index, up_limit);
375+
int hi_index = std::max<int>(k + current_width[k], lo_limit);
376+
hi_index = std::min<int>(k + current_width[k], up_limit);
403377
T_real temp = (background[lo_index] + background[hi_index]) / (T_real)2.0;
404378
if (background[k] > temp)
405379
{
@@ -418,6 +392,46 @@ DLL_EXPORT ArrayTr<T_real> snip_background(const Spectra<T_real> * const spectra
418392
}
419393

420394
// ----------------------------------------------------------------------------
395+
/*
396+
template<typename T_real>
397+
DLL_EXPORT ArrayTr<T_real> snip_background2(const Spectra<T_real> * const spectra, T_real width)
398+
{
399+
width = width * 100.0;
400+
assert(spectra != nullptr);
401+
ArrayTr<T_real> boxcar;
402+
boxcar.resize(5);
403+
boxcar.setConstant(5, 1.0);
404+
ArrayTr<T_real> background = convolve1d(*spectra, boxcar);
405+
background += 1.0;
406+
background = Eigen::sqrt(background);
407+
background = Eigen::log(Eigen::log(background + (T_real)1.0) + (T_real)1.0);
408+
409+
while (width >= 0.5)
410+
{
411+
for (int k = 0; k < background.size(); k++)
412+
{
413+
int lo = std::max<int>( k - (int)width, 0);
414+
int hi = std::min<int>( k + (int)width, background.size()-1);
415+
T_real temp = (background[lo] + background[hi]) / (T_real)2.0;
416+
if (background[k] > temp)
417+
{
418+
background[k] = temp;
419+
}
420+
}
421+
width = width / T_real(M_SQRT2); // window_rf
422+
}
423+
424+
background = Eigen::exp(Eigen::exp(background) - (T_real)1.0) - (T_real)1.0;
425+
background = Eigen::pow(background, 2.0);
426+
background -= 1.0;
427+
background = background.unaryExpr([](T_real v) { return std::isfinite(v) ? v : (T_real)0.0; });
428+
429+
return background;
430+
431+
}
432+
*/
433+
// ----------------------------------------------------------------------------
434+
421435

422436
template<typename T_real>
423437
using IO_Callback_Func_Def = std::function<void(size_t, size_t, size_t, size_t, size_t, Spectra<T_real>*, void*)>;

‎src/io/file/hdf5_io.h‎

Lines changed: 53 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -900,6 +900,59 @@ class DLL_EXPORT HDF5_IO
900900

901901
//-----------------------------------------------------------------------------
902902

903+
template<typename T_real>
904+
bool load_spectra_vol_apsu(std::string path, std::string &title, data_struct::Spectra_Volume<T_real>* spec_vol, data_struct::Scan_Info<T_real> &scan_info_edf, bool logerr = true)
905+
{
906+
907+
std::lock_guard<std::mutex> lock(_mutex);
908+
909+
std::stack<std::pair<hid_t, H5_OBJECTS> > close_map;
910+
hid_t file_id, xspres_grp_id, pos_grp_id, title_id;
911+
//hid_t root_grp_id, scanDim1_id, scanDim2_id, start_time_id;
912+
//H5G_info_t info;
913+
//char root_group_name[2048] = { 0 };
914+
915+
if (false == _open_h5_object(file_id, H5O_FILE, close_map, path, -1))
916+
return false;
917+
918+
//H5Gget_objname_by_idx(file_id, 0, &root_group_name[0], 2047);
919+
920+
//if (false == _open_h5_object(root_grp_id, H5O_GROUP, close_map, root_group_name, file_id))
921+
// return false;
922+
923+
if (false == _open_h5_object(pos_grp_id, H5O_GROUP, close_map, "/positions", file_id))
924+
return false;
925+
926+
if (false == _open_h5_object(xspres_grp_id, H5O_DATASET, close_map, "/xspress3", file_id))
927+
return false;
928+
/*
929+
if (false == _open_h5_object(scanDim2_id, H5O_DATASET, close_map, "scanDim_2", fluo_grp_id))
930+
return false;
931+
932+
if (false == _open_h5_object(title_id, H5O_DATASET, close_map, "title", root_grp_id))
933+
return false;
934+
935+
if (false == _open_h5_object(start_time_id, H5O_DATASET, close_map, "start_time", root_grp_id))
936+
return false;
937+
938+
char tmp_name[256] = { 0 };
939+
hid_t ftype = H5Dget_type(title_id);
940+
close_map.push({ ftype, H5O_DATATYPE });
941+
hid_t type = H5Tget_native_type(ftype, H5T_DIR_ASCEND);
942+
herr_t error = H5Dread(title_id, type, H5S_ALL, H5S_ALL, H5P_DEFAULT, (void*)&tmp_name[0]);
943+
944+
if (error == 0)
945+
{
946+
title = std::string(tmp_name);
947+
}
948+
*/
949+
950+
_close_h5_objects(close_map);
951+
return true;
952+
}
953+
954+
//-----------------------------------------------------------------------------
955+
903956
template<typename T_real>
904957
bool load_spectra_vol_esrf(std::string path, std::string &title, data_struct::Spectra_Volume<T_real>* spec_vol, data_struct::Scan_Info<T_real> &scan_info_edf, bool logerr = true)
905958
{

‎src/io/file/hl_file_io.h‎

Lines changed: 10 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -928,7 +928,16 @@ DLL_EXPORT bool load_spectra_volume(std::string dataset_directory,
928928
fullpath = dataset_directory + DIR_END_CHAR + dataset_file;
929929
std::string file_title;
930930
data_struct::Scan_Info<T_real> scan_info_edf;
931-
if(true == io::file::HDF5_IO::inst()->load_spectra_vol_esrf(fullpath, file_title, spectra_volume, scan_info_edf))
931+
// try new APS-U format
932+
if(true == io::file::HDF5_IO::inst()->load_spectra_vol_apsu(fullpath, file_title, spectra_volume, scan_info_edf))
933+
{
934+
// load scan rows from ./flyXRF
935+
936+
// load positions from ./positions
937+
938+
}
939+
// try ESRF dataset
940+
else if(true == io::file::HDF5_IO::inst()->load_spectra_vol_esrf(fullpath, file_title, spectra_volume, scan_info_edf))
932941
{
933942
std::string dset_folder = "";
934943
std::string base_name = dataset_file;

0 commit comments

Comments
 (0)