How to implement custom smoothing filter in DicomObjects
Smoothing can be subjective and we have enabled custom filtering in DicomObjects so users can define their own smoothing algorithms.
DicomObjects.NET version
public SharpDX.Vector4[] MakeSmoothingValues(out SharpDX.Vector2 offset, FilterMode mode, int size = 256)
{
SharpDX.Vector4[] data = new SharpDX.Vector4[size];
offset = new SharpDX.Vector2(-0.5F, -0.5F);
for (int i = 0; i < size; i++)
{
// Sample BSpline2 algorithm
float x = i / (float)(size - 1);
float x2 = x * x;
float x3 = x * x2;
if (mode == FilterMode.CustomValues)
{
data[i].W = (x3 - x2) / 2;
data[i].X = (-x3 + 2 * x2 - x) / 2;
data[i].Y = (3 * x3 - 5 * x2 + 2) / 2;
data[i].Z = (-3 * x3 + 4 * x2 + x) / 2;
}
else
{
offset = SharpDX.Vector2.Zero;
data[i] = new SharpDX.Vector4(0, 1, 0, 0);
}
}
return data;
}
dicomImage.MagnificationMode = FilterMode.CustomValues;
dicomImage.MinificationMode = FilterMode.CustomValues;
dicomImage.CustomSmoothingFunction = MakeSmoothingValues;
DicomObjects.COM version
Implement ICustomFilteringFunction interface
Due to COM limitation and lack of Delegate support, users must implment their own ICustomFilteringFunction interface.
In VB6 - create a new Class Module “MyFilter”:
Implements ICustomFilteringFunction
Public Function ICustomFilteringFunction_Filter(ByVal SmoothingSize As Long) As Variant
Dim size As Integer
Dim i As Integer
size = SmoothingSize * 4
Dim arr() As Single
ReDim arr(0 To size - 1)
Dim x, x2, x3 As Single
For i = 0 To SmoothingSize - 1
x = i / (SmoothingSize - 1)
x2 = x * x
x3 = x * x * x
arr(i * 4) = (x3 - x2) / 2 ' w
arr(i * 4 + 1) = (-x3 + 2 * x2 - x) / 2 ' x
arr(i * 4 + 2) = (3 * x3 - 5 * x2 + 2) / 2 ' y
arr(i * 4 + 3) = (-3 * x3 + 4 * x2 + x) / 2 ' z
Debug.Print "index:=" & i & " x:= " & (-x3 + 2 * x2 - x) / 2 & " y:= " & (3 * x3 - 5 * x2 + 2) / 2 & " z:= " & (-3 * x3 + 4 * x2 + x) / 2 & " w:= " & (x3 - x2) / 2
Next
ICustomFilteringFunction_Filter = arr ' returns SAFEARRAY(Single) in Variant
End Function
To set customer filter:
Dim f As New MyFilter
Viewer.CurrentImage.MagnificationMode = doFilterCustomValues
Viewer.CurrentImage.MinificationMode = doFilterCustomValues
Viewer.CurrentImage.SetCustomFilteringFunction f
In C++, create “FilterObject”
FilterObject.h
#pragma once
#include <windows.h>
#include <oleauto.h>
// IDispatch object exposing one late-bound method: "Filter(Long) As Variant"
// Returns Variant(SAFEARRAY of Single) of length SmoothingSize*4
class FilterObject final : public IDispatch
{
public:
// IUnknown
HRESULT STDMETHODCALLTYPE QueryInterface(REFIID riid, void** ppv) override;
ULONG STDMETHODCALLTYPE AddRef() override;
ULONG STDMETHODCALLTYPE Release() override;
// IDispatch
HRESULT STDMETHODCALLTYPE GetTypeInfoCount(UINT* pctinfo) override;
HRESULT STDMETHODCALLTYPE GetTypeInfo(UINT, LCID, ITypeInfo**) override;
HRESULT STDMETHODCALLTYPE GetIDsOfNames(REFIID, LPOLESTR* names, UINT cNames,
LCID, DISPID* dispids) override;
HRESULT STDMETHODCALLTYPE Invoke(DISPID dispid, REFIID, LCID,
WORD flags, DISPPARAMS* dp,
VARIANT* pRet, EXCEPINFO*,
UINT* pArgErr) override;
// Factory (refcount starts at 1)
static HRESULT Create(FilterObject** pp);
private:
FilterObject();
~FilterObject();
volatile LONG m_ref;
};
FilterObject.cpp
#include "stdafx.h"
#include "FilterObject.h"
#include <vector>
#include <sstream>
#include <iomanip>
#include <cwchar>
// --- helpers --------------------------------------------------------------
static void DebugPrintLine(long i, float x, float x2, float x3)
{
float X = (-x3 + 2.0f * x2 - x) / 2.0f;
float Y = (3.0f * x3 - 5.0f * x2 + 2.0f) / 2.0f;
float Z = (-3.0f * x3 + 4.0f * x2 + x) / 2.0f;
float W = (x3 - x2) / 2.0f;
std::wostringstream oss;
oss << L"index:= " << i
<< L" x:= " << std::fixed << std::setprecision(6) << X
<< L" y:= " << Y
<< L" z:= " << Z
<< L" w:= " << W << L"\r\n";
OutputDebugStringW(oss.str().c_str());
}
static HRESULT VectorToVariantR4(const std::vector<float>& v, VARIANT* pv)
{
if (!pv) return E_POINTER;
VariantInit(pv);
SAFEARRAYBOUND b{};
b.lLbound = 0;
b.cElements = static_cast<ULONG>(v.size());
SAFEARRAY* psa = SafeArrayCreate(VT_R4, 1, &b);
if (!psa) return E_OUTOFMEMORY;
float* pData = nullptr;
HRESULT hr = SafeArrayAccessData(psa, reinterpret_cast<void**>(&pData));
if (FAILED(hr)) { SafeArrayDestroy(psa); return hr; }
memcpy(pData, v.data(), v.size() * sizeof(float));
SafeArrayUnaccessData(psa);
pv->vt = VT_ARRAY | VT_R4;
pv->parray = psa;
return S_OK;
}
// --- FilterObject ---------------------------------------------------------
FilterObject::FilterObject() : m_ref(1) {}
FilterObject::~FilterObject() = default;
// IUnknown
HRESULT STDMETHODCALLTYPE FilterObject::QueryInterface(REFIID riid, void** ppv)
{
if (!ppv) return E_POINTER;
if (riid == IID_IUnknown || riid == IID_IDispatch) {
*ppv = static_cast<IDispatch*>(this);
AddRef();
return S_OK;
}
*ppv = nullptr;
return E_NOINTERFACE;
}
ULONG STDMETHODCALLTYPE FilterObject::AddRef()
{
return (ULONG)InterlockedIncrement(&m_ref);
}
ULONG STDMETHODCALLTYPE FilterObject::Release()
{
ULONG n = (ULONG)InterlockedDecrement(&m_ref);
if (n == 0) delete this;
return n;
}
// IDispatch
HRESULT STDMETHODCALLTYPE FilterObject::GetTypeInfoCount(UINT* pctinfo)
{
if (!pctinfo) return E_POINTER;
*pctinfo = 0;
return S_OK;
}
HRESULT STDMETHODCALLTYPE FilterObject::GetTypeInfo(UINT, LCID, ITypeInfo**)
{
return E_NOTIMPL;
}
HRESULT STDMETHODCALLTYPE FilterObject::GetIDsOfNames(REFIID, LPOLESTR* names, UINT cNames,
LCID, DISPID* dispids)
{
if (!names || !dispids) return E_POINTER;
if (cNames < 1) return E_INVALIDARG;
if (_wcsicmp(names[0], L"ICustomFilteringFunction_Filter") == 0) {
dispids[0] = 1; // our chosen DISPID
return S_OK;
}
return DISP_E_UNKNOWNNAME;
}
HRESULT STDMETHODCALLTYPE FilterObject::Invoke(DISPID dispid, REFIID, LCID,
WORD flags, DISPPARAMS* dp,
VARIANT* pRet, EXCEPINFO*,
UINT* pArgErr)
{
if (dispid != 1 || !(flags & DISPATCH_METHOD))
return DISP_E_MEMBERNOTFOUND;
if (!dp || dp->cArgs != 1) {
if (pArgErr) *pArgErr = 0;
return DISP_E_BADPARAMCOUNT;
}
// Coerce arg0 to LONG (SmoothingSize)
VARIANT vIn; VariantInit(&vIn);
HRESULT hr = VariantChangeType(&vIn, &dp->rgvarg[0], 0, VT_I4);
if (FAILED(hr)) {
if (pArgErr) *pArgErr = 0;
return DISP_E_TYPEMISMATCH;
}
LONG smoothingSize = vIn.lVal;
if (smoothingSize < 0) smoothingSize = 0;
std::vector<float> out(static_cast<size_t>(smoothingSize) * 4u);
const bool movingAverage = false;
for (LONG i = 0; i < smoothingSize; ++i) {
size_t base = static_cast<size_t>(i) * 4u;
if (movingAverage) {
out[base + 0] = 0.25f;
out[base + 1] = 0.25f;
out[base + 2] = 0.25f;
out[base + 3] = 0.25f;
float denom = (smoothingSize > 1) ? float(smoothingSize - 1) : 1.0f;
float x = (smoothingSize > 1) ? (float(i) / denom) : 0.0f;
float x2 = x * x;
float x3 = x2 * x;
DebugPrintLine(i, x, x2, x3);
}
else { // sample BSpline2
float denom = (smoothingSize > 1) ? float(smoothingSize - 1) : 1.0f;
float x = (smoothingSize > 1) ? (float(i) / denom) : 0.0f;
float x2 = x * x;
float x3 = x2 * x;
float w = (x3 - x2) / 2.0f;
float X = (-x3 + 2.0f * x2 - x) / 2.0f;
float Y = (3.0f * x3 - 5.0f * x2 + 2.0f) / 2.0f;
float Z = (-3.0f * x3 + 4.0f * x2 + x) / 2.0f;
out[base + 0] = w;
out[base + 1] = X;
out[base + 2] = Y;
out[base + 3] = Z;
DebugPrintLine(i, x, x2, x3);
}
}
if (pRet) return VectorToVariantR4(out, pRet);
return S_OK;
}
// Factory
HRESULT FilterObject::Create(FilterObject** pp)
{
if (!pp) return E_POINTER;
*pp = new (std::nothrow) FilterObject();
return (*pp) ? S_OK : E_OUTOFMEMORY;
}
FilterObject* pFilter = nullptr;
if (SUCCEEDED(FilterObject::Create(&pFilter)))
{
ActiveImage->SetCustomFilteringFunction(static_cast<IDispatch*>(pFilter));
ActiveImage->UnsharpEnhancement = -0.5f;
ActiveImage->UnsharpLength = 2;
pFilter->Release();
}