Skip to content

Commit

Permalink
Merge pull request #195 from aglowacki/master
Browse files Browse the repository at this point in the history 
Refactor snip
  • Loading branch information
@aglowacki
aglowacki authored Jul 15, 2024
2 parents a457ec6 + d71a5c7 commit bb22360
Show file tree
Hide file tree
Showing 3 changed files with 116 additions and 40 deletions.
92 changes: 53 additions & 39 deletions src/data_struct/spectra.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,7 @@ POSSIBILITY OF SUCH DAMAGE.
#include "core/defines.h"
#include <Eigen/Core>
#include <vector>
#include <algorithm>
#include <functional>

namespace data_struct
Expand Down Expand Up @@ -254,6 +255,7 @@ class Spectra : public ArrayTr<_T>
TEMPLATE_CLASS_DLL_EXPORT Spectra<float>;
TEMPLATE_CLASS_DLL_EXPORT Spectra<double>;

//-----------------------------------------------------------------------------
// ----------------------------------------------------------------------------

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

// FIRST SNIPPING
int no_iterations = 2;

int max_of_xmin = (std::max)(xmin, (T_real)0.0);
int min_of_xmax = (std::min)(xmax, T_real(spectra->size() - 1));
int lo_limit = std::max<int>(xmin, 0);
int up_limit = std::min<int>(xmax, (spectra->size() - 1));
for (int j = 0; j < no_iterations; j++)
{
for (long int k = 0; k < background.size(); k++)
{
long int lo_index = k - current_width[k];
long int hi_index = k + current_width[k];
if (lo_index < max_of_xmin)
{
lo_index = max_of_xmin;
}
if (lo_index > min_of_xmax)
{
lo_index = min_of_xmax;
}
if (hi_index > min_of_xmax)
{
hi_index = min_of_xmax;
}
if (hi_index < max_of_xmin)
{
hi_index = max_of_xmin;
}
int lo_index = std::max<int>(k - current_width[k], lo_limit);
lo_index = std::min<int>(lo_index, up_limit);
int hi_index = std::max<int>(k + current_width[k], lo_limit);
hi_index = std::min<int>(k + current_width[k], up_limit);
T_real temp = (background[lo_index] + background[hi_index]) / (T_real)2.0;
if (background[k] > temp)
{
Expand All @@ -382,24 +370,10 @@ DLL_EXPORT ArrayTr<T_real> snip_background(const Spectra<T_real> * const spectra
{
for (long int k = 0; k < background.size(); k++)
{
long int lo_index = k - current_width[k];
long int hi_index = k + current_width[k];
if (lo_index < max_of_xmin)
{
lo_index = max_of_xmin;
}
if (lo_index > min_of_xmax)
{
lo_index = min_of_xmax;
}
if (hi_index > min_of_xmax)
{
hi_index = min_of_xmax;
}
if (hi_index < max_of_xmin)
{
hi_index = max_of_xmin;
}
int lo_index = std::max<int>(k - current_width[k], lo_limit);
lo_index = std::min<int>(lo_index, up_limit);
int hi_index = std::max<int>(k + current_width[k], lo_limit);
hi_index = std::min<int>(k + current_width[k], up_limit);
T_real temp = (background[lo_index] + background[hi_index]) / (T_real)2.0;
if (background[k] > temp)
{
Expand All @@ -418,6 +392,46 @@ DLL_EXPORT ArrayTr<T_real> snip_background(const Spectra<T_real> * const spectra
}

// ----------------------------------------------------------------------------
/*
template<typename T_real>
DLL_EXPORT ArrayTr<T_real> snip_background2(const Spectra<T_real> * const spectra, T_real width)
{
width = width * 100.0;
assert(spectra != nullptr);
ArrayTr<T_real> boxcar;
boxcar.resize(5);
boxcar.setConstant(5, 1.0);
ArrayTr<T_real> background = convolve1d(*spectra, boxcar);
background += 1.0;
background = Eigen::sqrt(background);
background = Eigen::log(Eigen::log(background + (T_real)1.0) + (T_real)1.0);
while (width >= 0.5)
{
for (int k = 0; k < background.size(); k++)
{
int lo = std::max<int>( k - (int)width, 0);
int hi = std::min<int>( k + (int)width, background.size()-1);
T_real temp = (background[lo] + background[hi]) / (T_real)2.0;
if (background[k] > temp)
{
background[k] = temp;
}
}
width = width / T_real(M_SQRT2); // window_rf
}
background = Eigen::exp(Eigen::exp(background) - (T_real)1.0) - (T_real)1.0;
background = Eigen::pow(background, 2.0);
background -= 1.0;
background = background.unaryExpr([](T_real v) { return std::isfinite(v) ? v : (T_real)0.0; });
return background;
}
*/
// ----------------------------------------------------------------------------


template<typename T_real>
using IO_Callback_Func_Def = std::function<void(size_t, size_t, size_t, size_t, size_t, Spectra<T_real>*, void*)>;
Expand Down
53 changes: 53 additions & 0 deletions src/io/file/hdf5_io.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -900,6 +900,59 @@ class DLL_EXPORT HDF5_IO

//-----------------------------------------------------------------------------

template<typename T_real>
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)
{

std::lock_guard<std::mutex> lock(_mutex);

std::stack<std::pair<hid_t, H5_OBJECTS> > close_map;
hid_t file_id, xspres_grp_id, pos_grp_id, title_id;
//hid_t root_grp_id, scanDim1_id, scanDim2_id, start_time_id;
//H5G_info_t info;
//char root_group_name[2048] = { 0 };

if (false == _open_h5_object(file_id, H5O_FILE, close_map, path, -1))
return false;

//H5Gget_objname_by_idx(file_id, 0, &root_group_name[0], 2047);

//if (false == _open_h5_object(root_grp_id, H5O_GROUP, close_map, root_group_name, file_id))
// return false;

if (false == _open_h5_object(pos_grp_id, H5O_GROUP, close_map, "/positions", file_id))
return false;

if (false == _open_h5_object(xspres_grp_id, H5O_DATASET, close_map, "/xspress3", file_id))
return false;
/*
if (false == _open_h5_object(scanDim2_id, H5O_DATASET, close_map, "scanDim_2", fluo_grp_id))
return false;
if (false == _open_h5_object(title_id, H5O_DATASET, close_map, "title", root_grp_id))
return false;
if (false == _open_h5_object(start_time_id, H5O_DATASET, close_map, "start_time", root_grp_id))
return false;
char tmp_name[256] = { 0 };
hid_t ftype = H5Dget_type(title_id);
close_map.push({ ftype, H5O_DATATYPE });
hid_t type = H5Tget_native_type(ftype, H5T_DIR_ASCEND);
herr_t error = H5Dread(title_id, type, H5S_ALL, H5S_ALL, H5P_DEFAULT, (void*)&tmp_name[0]);
if (error == 0)
{
title = std::string(tmp_name);
}
*/

_close_h5_objects(close_map);
return true;
}

//-----------------------------------------------------------------------------

template<typename T_real>
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)
{
Expand Down
11 changes: 10 additions & 1 deletion src/io/file/hl_file_io.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -928,7 +928,16 @@ DLL_EXPORT bool load_spectra_volume(std::string dataset_directory,
fullpath = dataset_directory + DIR_END_CHAR + dataset_file;
std::string file_title;
data_struct::Scan_Info<T_real> scan_info_edf;
if(true == io::file::HDF5_IO::inst()->load_spectra_vol_esrf(fullpath, file_title, spectra_volume, scan_info_edf))
// try new APS-U format
if(true == io::file::HDF5_IO::inst()->load_spectra_vol_apsu(fullpath, file_title, spectra_volume, scan_info_edf))
{
// load scan rows from ./flyXRF

// load positions from ./positions

}
// try ESRF dataset
else if(true == io::file::HDF5_IO::inst()->load_spectra_vol_esrf(fullpath, file_title, spectra_volume, scan_info_edf))
{
std::string dset_folder = "";
std::string base_name = dataset_file;
Expand Down

0 comments on commit bb22360

Please sign in to comment.