Mercurial > hg > orthanc
diff OrthancServer/Sources/SliceOrdering.cpp @ 4044:d25f4c0fa160 framework
splitting code into OrthancFramework and OrthancServer
author | Sebastien Jodogne <s.jodogne@gmail.com> |
---|---|
date | Wed, 10 Jun 2020 20:30:34 +0200 |
parents | OrthancServer/SliceOrdering.cpp@104e27133ebd |
children | 05b8fd21089c |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/OrthancServer/Sources/SliceOrdering.cpp Wed Jun 10 20:30:34 2020 +0200 @@ -0,0 +1,491 @@ +/** + * Orthanc - A Lightweight, RESTful DICOM Store + * Copyright (C) 2012-2016 Sebastien Jodogne, Medical Physics + * Department, University Hospital of Liege, Belgium + * Copyright (C) 2017-2020 Osimis S.A., Belgium + * + * This program is free software: you can redistribute it and/or + * modify it under the terms of the GNU General Public License as + * published by the Free Software Foundation, either version 3 of the + * License, or (at your option) any later version. + * + * In addition, as a special exception, the copyright holders of this + * program give permission to link the code of its release with the + * OpenSSL project's "OpenSSL" library (or with modified versions of it + * that use the same license as the "OpenSSL" library), and distribute + * the linked executables. You must obey the GNU General Public License + * in all respects for all of the code used other than "OpenSSL". If you + * modify file(s) with this exception, you may extend this exception to + * your version of the file(s), but you are not obligated to do so. If + * you do not wish to do so, delete this exception statement from your + * version. If you delete this exception statement from all source files + * in the program, then also delete it here. + * + * This program is distributed in the hope that it will be useful, but + * WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + **/ + + +#include "PrecompiledHeadersServer.h" +#include "SliceOrdering.h" + +#include "../Core/Logging.h" +#include "../Core/Toolbox.h" +#include "ServerEnumerations.h" +#include "ServerIndex.h" + +#include <algorithm> +#include <boost/lexical_cast.hpp> +#include <boost/noncopyable.hpp> + + +namespace Orthanc +{ + static bool TokenizeVector(std::vector<float>& result, + const std::string& value, + unsigned int expectedSize) + { + std::vector<std::string> tokens; + Toolbox::TokenizeString(tokens, value, '\\'); + + if (tokens.size() != expectedSize) + { + return false; + } + + result.resize(tokens.size()); + + for (size_t i = 0; i < tokens.size(); i++) + { + try + { + result[i] = boost::lexical_cast<float>(tokens[i]); + } + catch (boost::bad_lexical_cast&) + { + return false; + } + } + + return true; + } + + + static bool TokenizeVector(std::vector<float>& result, + const DicomMap& map, + const DicomTag& tag, + unsigned int expectedSize) + { + const DicomValue* value = map.TestAndGetValue(tag); + + if (value == NULL || + value->IsNull() || + value->IsBinary()) + { + return false; + } + else + { + return TokenizeVector(result, value->GetContent(), expectedSize); + } + } + + + static bool IsCloseToZero(double x) + { + return fabs(x) < 10.0 * std::numeric_limits<float>::epsilon(); + } + + + bool SliceOrdering::ComputeNormal(Vector& normal, + const DicomMap& dicom) + { + std::vector<float> cosines; + + if (TokenizeVector(cosines, dicom, DICOM_TAG_IMAGE_ORIENTATION_PATIENT, 6)) + { + assert(cosines.size() == 6); + normal[0] = cosines[1] * cosines[5] - cosines[2] * cosines[4]; + normal[1] = cosines[2] * cosines[3] - cosines[0] * cosines[5]; + normal[2] = cosines[0] * cosines[4] - cosines[1] * cosines[3]; + return true; + } + else + { + return false; + } + } + + + bool SliceOrdering::IsParallelOrOpposite(const Vector& u, + const Vector& v) + { + // Check out "GeometryToolbox::IsParallelOrOpposite()" in Stone of + // Orthanc for explanations + const double u1 = u[0]; + const double u2 = u[1]; + const double u3 = u[2]; + const double normU = sqrt(u1 * u1 + u2 * u2 + u3 * u3); + + const double v1 = v[0]; + const double v2 = v[1]; + const double v3 = v[2]; + const double normV = sqrt(v1 * v1 + v2 * v2 + v3 * v3); + + if (IsCloseToZero(normU * normV)) + { + return false; + } + else + { + const double cosAngle = (u1 * v1 + u2 * v2 + u3 * v3) / (normU * normV); + + return (IsCloseToZero(cosAngle - 1.0) || // Close to +1: Parallel, non-opposite + IsCloseToZero(fabs(cosAngle) - 1.0)); // Close to -1: Parallel, opposite + } + } + + + struct SliceOrdering::Instance : public boost::noncopyable + { + private: + std::string instanceId_; + bool hasPosition_; + Vector position_; + bool hasNormal_; + Vector normal_; + bool hasIndexInSeries_; + size_t indexInSeries_; + unsigned int framesCount_; + + public: + Instance(ServerIndex& index, + const std::string& instanceId) : + instanceId_(instanceId), + framesCount_(1) + { + DicomMap instance; + if (!index.GetMainDicomTags(instance, instanceId, ResourceType_Instance, ResourceType_Instance)) + { + throw OrthancException(ErrorCode_UnknownResource); + } + + const DicomValue* frames = instance.TestAndGetValue(DICOM_TAG_NUMBER_OF_FRAMES); + if (frames != NULL && + !frames->IsNull() && + !frames->IsBinary()) + { + try + { + framesCount_ = boost::lexical_cast<unsigned int>(frames->GetContent()); + } + catch (boost::bad_lexical_cast&) + { + } + } + + std::vector<float> tmp; + hasPosition_ = TokenizeVector(tmp, instance, DICOM_TAG_IMAGE_POSITION_PATIENT, 3); + + if (hasPosition_) + { + position_[0] = tmp[0]; + position_[1] = tmp[1]; + position_[2] = tmp[2]; + } + + hasNormal_ = ComputeNormal(normal_, instance); + + std::string s; + hasIndexInSeries_ = false; + + try + { + if (index.LookupMetadata(s, instanceId, MetadataType_Instance_IndexInSeries)) + { + indexInSeries_ = boost::lexical_cast<size_t>(s); + hasIndexInSeries_ = true; + } + } + catch (boost::bad_lexical_cast&) + { + } + } + + const std::string& GetIdentifier() const + { + return instanceId_; + } + + bool HasPosition() const + { + return hasPosition_; + } + + float ComputeRelativePosition(const Vector& normal) const + { + assert(HasPosition()); + return (normal[0] * position_[0] + + normal[1] * position_[1] + + normal[2] * position_[2]); + } + + bool HasIndexInSeries() const + { + return hasIndexInSeries_; + } + + size_t GetIndexInSeries() const + { + assert(HasIndexInSeries()); + return indexInSeries_; + } + + unsigned int GetFramesCount() const + { + return framesCount_; + } + + bool HasNormal() const + { + return hasNormal_; + } + + const Vector& GetNormal() const + { + assert(hasNormal_); + return normal_; + } + }; + + + class SliceOrdering::PositionComparator + { + private: + const Vector& normal_; + + public: + PositionComparator(const Vector& normal) : normal_(normal) + { + } + + int operator() (const Instance* a, + const Instance* b) const + { + return a->ComputeRelativePosition(normal_) < b->ComputeRelativePosition(normal_); + } + }; + + + bool SliceOrdering::IndexInSeriesComparator(const SliceOrdering::Instance* a, + const SliceOrdering::Instance* b) + { + return a->GetIndexInSeries() < b->GetIndexInSeries(); + } + + + void SliceOrdering::ComputeNormal() + { + DicomMap series; + if (!index_.GetMainDicomTags(series, seriesId_, ResourceType_Series, ResourceType_Series)) + { + throw OrthancException(ErrorCode_UnknownResource); + } + + hasNormal_ = ComputeNormal(normal_, series); + } + + + void SliceOrdering::CreateInstances() + { + std::list<std::string> instancesId; + index_.GetChildren(instancesId, seriesId_); + + instances_.reserve(instancesId.size()); + for (std::list<std::string>::const_iterator + it = instancesId.begin(); it != instancesId.end(); ++it) + { + instances_.push_back(new Instance(index_, *it)); + } + } + + + bool SliceOrdering::SortUsingPositions() + { + if (instances_.size() <= 1) + { + // One single instance: It is sorted by default + return true; + } + + if (!hasNormal_) + { + return false; + } + + for (size_t i = 0; i < instances_.size(); i++) + { + assert(instances_[i] != NULL); + + if (!instances_[i]->HasPosition() || + (instances_[i]->HasNormal() && + !IsParallelOrOpposite(instances_[i]->GetNormal(), normal_))) + { + return false; + } + } + + PositionComparator comparator(normal_); + std::sort(instances_.begin(), instances_.end(), comparator); + + float a = instances_[0]->ComputeRelativePosition(normal_); + for (size_t i = 1; i < instances_.size(); i++) + { + float b = instances_[i]->ComputeRelativePosition(normal_); + + if (std::fabs(b - a) <= 10.0f * std::numeric_limits<float>::epsilon()) + { + // Not enough space between two slices along the normal of the volume + return false; + } + + a = b; + } + + // This is a 3D volume + isVolume_ = true; + return true; + } + + + bool SliceOrdering::SortUsingIndexInSeries() + { + if (instances_.size() <= 1) + { + // One single instance: It is sorted by default + return true; + } + + for (size_t i = 0; i < instances_.size(); i++) + { + assert(instances_[i] != NULL); + if (!instances_[i]->HasIndexInSeries()) + { + return false; + } + } + + std::sort(instances_.begin(), instances_.end(), IndexInSeriesComparator); + + for (size_t i = 1; i < instances_.size(); i++) + { + if (instances_[i - 1]->GetIndexInSeries() == instances_[i]->GetIndexInSeries()) + { + // The current "IndexInSeries" occurs 2 times: Not a proper ordering + LOG(WARNING) << "This series contains 2 slices with the same index, trying to display it anyway"; + break; + } + } + + return true; + } + + + SliceOrdering::SliceOrdering(ServerIndex& index, + const std::string& seriesId) : + index_(index), + seriesId_(seriesId), + isVolume_(false) + { + ComputeNormal(); + CreateInstances(); + + if (!SortUsingPositions() && + !SortUsingIndexInSeries()) + { + throw OrthancException(ErrorCode_CannotOrderSlices, + "Unable to order the slices of series " + seriesId); + } + } + + + SliceOrdering::~SliceOrdering() + { + for (std::vector<Instance*>::iterator + it = instances_.begin(); it != instances_.end(); ++it) + { + if (*it != NULL) + { + delete *it; + } + } + } + + + const std::string& SliceOrdering::GetInstanceId(size_t index) const + { + if (index >= instances_.size()) + { + throw OrthancException(ErrorCode_ParameterOutOfRange); + } + else + { + return instances_[index]->GetIdentifier(); + } + } + + + unsigned int SliceOrdering::GetFramesCount(size_t index) const + { + if (index >= instances_.size()) + { + throw OrthancException(ErrorCode_ParameterOutOfRange); + } + else + { + return instances_[index]->GetFramesCount(); + } + } + + + void SliceOrdering::Format(Json::Value& result) const + { + result = Json::objectValue; + result["Type"] = (isVolume_ ? "Volume" : "Sequence"); + + Json::Value tmp = Json::arrayValue; + for (size_t i = 0; i < GetInstancesCount(); i++) + { + tmp.append(GetBasePath(ResourceType_Instance, GetInstanceId(i)) + "/file"); + } + + result["Dicom"] = tmp; + + Json::Value slicesShort = Json::arrayValue; + + tmp.clear(); + for (size_t i = 0; i < GetInstancesCount(); i++) + { + std::string base = GetBasePath(ResourceType_Instance, GetInstanceId(i)); + for (size_t j = 0; j < GetFramesCount(i); j++) + { + tmp.append(base + "/frames/" + boost::lexical_cast<std::string>(j)); + } + + Json::Value tmp2 = Json::arrayValue; + tmp2.append(GetInstanceId(i)); + tmp2.append(0); + tmp2.append(GetFramesCount(i)); + + slicesShort.append(tmp2); + } + + result["Slices"] = tmp; + result["SlicesShort"] = slicesShort; + } +}