elemental_profile.cpp

Go to the documentation of this file.

 
#include "elemental_profile.h"
#include <iostream>
#include "Utilities.h"
#include <QStringList>
#include <QFile>
 
 
// =============================================================================
// CONSTRUCTION AND ASSIGNMENT
// =============================================================================
 
Elemental_Profile::Elemental_Profile()
    : map<string, double>(),
    Interface(),
    included_in_analysis_(true),
    marked_elements_()
{
    // All initialization done in initializer list
    // Empty body is fine - demonstrates proper C++ style
}
 
 
Elemental_Profile::Elemental_Profile(const Elemental_Profile& other)
    : map<string, double>(other),  // Copy base class (element concentrations)
    Interface(other),              // Copy Interface base class
    included_in_analysis_(other.included_in_analysis_),
    marked_elements_(other.marked_elements_)
{
    // All copying done in initializer list - efficient and exception-safe
}
 
 
/*/
*@brief Copy assignment operator
*
* Assigns the contents of another profile to this one, replacing all
* existing data.Handles self - assignment safely.
*
*@param other Profile to assign from
* @return Reference to * this for chaining(e.g., a = b = c)
*
* @post* this contains copy of other's data
* @post Previous contents of * this are replaced
* @post Self - assignment(a = a) is safe
*
*@note Assignment order is important : Interface first, then members,
* then base class to ensure proper cleanup if exceptions occur
*/
Elemental_Profile & Elemental_Profile::operator=(const Elemental_Profile & other)
{
    // Check for self-assignment (important for performance and correctness)
    if (this == &other) {
        return *this;
    }
 
    // Assign in proper order:
    // 1. Interface base class
    Interface::operator=(other);
 
    // 2. Member variables
    included_in_analysis_ = other.included_in_analysis_;
    marked_elements_ = other.marked_elements_;
 
    // 3. map base class (element concentrations)
    //    Done last because it might throw, and we want other data copied first
    map<string, double>::operator=(other);
 
    return *this;
}
 
Elemental_Profile Elemental_Profile::CreateAnalysisProfile(
    const map<string, element_information>* element_info) const
{
    // If no metadata provided, return complete copy
    if (element_info == nullptr) {
        return *this;
    }
 
    // Create filtered profile
    Elemental_Profile filtered_profile;
 
    // Copy analysis inclusion flag
    filtered_profile.included_in_analysis_ = this->included_in_analysis_;
 
    // Iterate through all elements in current profile
    for (const auto& element_entry : *this) {
        const string& element_name = element_entry.first;
        const double concentration = element_entry.second;
 
        // Check if element exists in metadata
        auto info_iter = element_info->find(element_name);
        if (info_iter == element_info->end()) {
            // Element not in metadata - skip it
            continue;
        }
 
        const element_information& info = info_iter->second;
 
        // Apply inclusion criteria
        bool should_include = info.include_in_analysis &&
            info.Role != element_information::role::do_not_include &&
            info.Role != element_information::role::particle_size &&
            info.Role != element_information::role::organic_carbon;
 
        if (should_include) {
            // Add element to filtered profile
            filtered_profile[element_name] = concentration;
 
            // Preserve marked status if element was marked
            auto marked_iter = marked_elements_.find(element_name);
            if (marked_iter != marked_elements_.end() && marked_iter->second) {
                filtered_profile.marked_elements_[element_name] = true;
            }
        }
    }
 
    return filtered_profile;
}
 
Elemental_Profile Elemental_Profile::CreateCorrectedProfile(
    bool exclude_elements,
    bool apply_corrections,
    const vector<double>& reference_om_size,
    const MultipleLinearRegressionSet* regression_models,
    const map<string, element_information>* element_info) const
{
    // Start with current profile
    Elemental_Profile result = *this;
 
    // Step 1: Apply OM/size corrections if requested
    if (apply_corrections) {
        if (regression_models == nullptr) {
            // Log warning or throw exception in production code
            std::cerr << "Warning: apply_corrections=true but regression_models is null. "
                << "Skipping corrections." << std::endl;
        }
        else {
            result = result.ApplyOrganicMatterAndSizeCorrections(
                reference_om_size,
                regression_models,
                element_info
            );
        }
    }
 
    // Step 2: Filter elements if requested
    if (exclude_elements) {
        if (element_info == nullptr) {
            // Log warning or throw exception in production code
            std::cerr << "Warning: exclude_elements=true but element_info is null. "
                << "Skipping filtering." << std::endl;
        }
        else {
            result = result.CreateAnalysisProfile(element_info);
        }
    }
 
    return result;
}
 
Elemental_Profile Elemental_Profile::ExtractChemicalElements(
    const map<string, element_information>* element_info,
    bool include_isotopes) const
{
    // Validate input - element_info is required
    if (element_info == nullptr) {
        // For now, log error and return empty profile
        // Future: throw std::invalid_argument("element_info cannot be null")
        std::cerr << "Error: ExtractChemicalElements called with null element_info. "
            << "Returning empty profile." << std::endl;
        return Elemental_Profile();
    }
 
    Elemental_Profile extracted_profile;
 
    // Iterate through all elements in current profile
    for (const auto& [element_name, concentration] : *this) {
 
        // Check if element exists in metadata
        auto info_iter = element_info->find(element_name);
        if (info_iter == element_info->end()) {
            // Element not in metadata - skip it
            // Could log warning: std::cerr << "Warning: Element '" << element_name 
            //                              << "' not found in metadata. Skipping.\n";
            continue;
        }
 
        const element_information& info = info_iter->second;
 
        // Check if element should be included based on role
        bool should_include = false;
 
        if (info.Role == element_information::role::element) {
            // Always include chemical elements
            should_include = true;
        }
        else if (info.Role == element_information::role::isotope) {
            // Include isotopes only if requested
            should_include = include_isotopes;
        }
        // All other roles (organic_carbon, particle_size, do_not_include) are excluded
 
        if (should_include) {
            extracted_profile[element_name] = concentration;
        }
    }
 
    return extracted_profile;
}
 
Elemental_Profile Elemental_Profile::ExtractElements(
    const vector<string>& element_names) const
{
    Elemental_Profile extracted_profile;
 
    // Reserve space if we expect most elements to exist (optimization)
    // Note: Can't reserve in map, but this shows intent
 
    // Iterate through requested elements
    for (const string& element_name : element_names) {
        // Check if element exists in current profile
        auto it = find(element_name);
        if (it != end()) {
            // Element exists - add to extracted profile
            extracted_profile[element_name] = it->second;
        }
        // Element doesn't exist - silently skip
    }
 
    return extracted_profile;
}
 
double Elemental_Profile::GetValue(const string& element_name) const
{
    auto it = find(element_name);
    if (it == end()) {
        throw std::out_of_range(
            "Element '" + element_name + "' does not exist in profile"
        );
    }
    return it->second;
}
 
 
 bool Elemental_Profile::SetValue(const string& element_name, double value)
{
    auto it = find(element_name);
    if (it == end()) {
        return false;
    }
 
    it->second = value;
    return true;
}
 
 bool Elemental_Profile::SetMarked(const string& element_name, bool marked)
 {
     if (find(element_name) == end()) {
         return false;
     }
 
     marked_elements_[element_name] = marked;
     return true;
 }
 
 bool Elemental_Profile::IsMarked(const string& element_name) const
 {
     if (find(element_name) == end()) {
         return false;
     }
 
     auto it = marked_elements_.find(element_name);
     return (it != marked_elements_.end() && it->second);
 }
 
 
 bool Elemental_Profile::AppendElement(const string& element_name, double value)
 {
     if (find(element_name) != end()) {
         return false;
     }
 
     (*this)[element_name] = value;
     marked_elements_[element_name] = false;
     return true;
 }
 
 vector<double> Elemental_Profile::GetAllValues() const
 {
     vector<double> values;
     values.reserve(size());
 
     for (const auto& [element_name, concentration] : *this) {
         values.push_back(concentration);
     }
 
     return values;
 }
 
 string Elemental_Profile::ToString() const
 {
     string output;
 
     for (const auto& [element_name, concentration] : *this) {
         output += element_name + ":" + aquiutils::numbertostring(concentration) + "\n";
     }
 
     return output;
 }
 
 QTableWidget* Elemental_Profile::ToTable()
 {
     QTableWidget* table = new QTableWidget();
     table->setEditTriggers(QAbstractItemView::NoEditTriggers);
     table->setColumnCount(1);
     table->setRowCount(size());
 
     QStringList headers;
     QStringList element_names;
 
     int row = 0;
     for (const auto& [element_name, concentration] : *this) {
         element_names << QString::fromStdString(element_name);
         table->setItem(row, 0, new QTableWidgetItem(QString::number(concentration)));
 
         // Highlight values outside limits if enabled
         if (highlightoutsideoflimit &&
             (concentration > highlimit || concentration < lowlimit)) {
             table->item(row, 0)->setForeground(QColor(Qt::red));
         }
 
         row++;
     }
 
     headers << YAxisLabel();
     table->setHorizontalHeaderLabels(headers);
     table->setVerticalHeaderLabels(element_names);
 
     return table;
 }
 
 QJsonObject Elemental_Profile::toJsonObject() const
 {
     QJsonObject json_object;
     QJsonObject element_contents;
 
     json_object["Include"] = IsIncludedInAnalysis();
 
     for (const auto& [element_name, concentration] : *this) {
         element_contents[QString::fromStdString(element_name)] = concentration;
     }
 
     json_object["Element Contents"] = element_contents;
     json_object["XAxisLabel"] = XAxisLabel();
     json_object["YAxisLabel"] = YAxisLabel();
 
     return json_object;
 }
 
 bool Elemental_Profile::ReadFromJsonObject(const QJsonObject& json_object)
 {
     clear();
 
     if (json_object.contains("Include")) {
         // Full format with metadata
         SetIncludedInAnalysis(json_object["Include"].toBool());
 
         QJsonObject contents = json_object["Element Contents"].toObject();
         for (const QString& key : contents.keys()) {
             (*this)[key.toStdString()] = contents[key].toDouble();
         }
     }
     else {
         // Legacy format - elements directly in root
         for (const QString& key : json_object.keys()) {
             (*this)[key.toStdString()] = json_object[key].toDouble();
         }
     }
 
     if (json_object.contains("XAxisLabel")) {
         SetXAxisLabel(json_object["XAxisLabel"].toString());
     }
 
     if (json_object.contains("YAxisLabel")) {
         SetYAxisLabel(json_object["YAxisLabel"].toString());  // Fixed: was SetXAxisLabel
     }
 
     return true;
 }
 
 
 bool Elemental_Profile::Read(const QStringList& string_list)
 {
     clear();
 
     for (const QString& line : string_list) {
         QStringList parts = line.split(":");
         if (parts.size() > 1) {
             AppendElement(parts[0].toStdString(), parts[1].toDouble());
         }
     }
 
     return true;
 }
 
 bool Elemental_Profile::writetofile(QFile* file)
 {
     file->write(QString::fromStdString(ToString()).toUtf8());
     return true;
 }
 
 double Elemental_Profile::GetMaximum() const
 {
     if (empty()) {
         return -1e12;
     }
 
     double max_value = -1e12;
     for (const auto& [element_name, concentration] : *this) {
         if (concentration > max_value) {
             max_value = concentration;
         }
     }
 
     return max_value;
 }
 
 double Elemental_Profile::GetMinimum() const
 {
     if (empty()) {
         return 1e12;
     }
 
     double min_value = 1e12;
     for (const auto& [element_name, concentration] : *this) {
         if (concentration < min_value) {
             min_value = concentration;
         }
     }
 
     return min_value;
 }
 
 vector<string> Elemental_Profile::GetElementNames() const
 {
     vector<string> names;
     names.reserve(size());
 
     for (const auto& [element_name, concentration] : *this) {
         names.push_back(element_name);
     }
 
     return names;
 }
 
 bool Elemental_Profile::Contains(const string& element_name) const
 {
     return (find(element_name) != end());
 }
 
 
 Elemental_Profile Elemental_Profile::ApplyOrganicMatterAndSizeCorrections(
     const vector<double>& reference_om_size,
     const MultipleLinearRegressionSet* regression_models,
     const map<string, element_information>* element_info) const
 {
     Elemental_Profile corrected_profile;
 
     for (const auto& [element_name, concentration] : *this) {
         const element_information& info = element_info->at(element_name);
 
         // Skip organic carbon and particle size indicators
         if (info.Role == element_information::role::organic_carbon ||
             info.Role == element_information::role::particle_size) {
             continue;
         }
 
         // Start with original concentration
         corrected_profile[element_name] = concentration;
 
         // Apply corrections if regression model exists for this element
         if (regression_models->count(element_name) == 0) {
             continue;
         }
 
         const MultipleLinearRegression* mlr = &regression_models->at(element_name);
         const auto& var_names = mlr->GetIndependentVariableNames();
         const auto& coefficients = mlr->CoefficientsIntercept();
         bool is_linear = (mlr->Equation() == regression_form::linear);
 
         // Apply first correction (typically organic matter)
         if (mlr->Effective(0) && var_names[0] != element_name) {
             double reference = reference_om_size[0];
             double measured = this->at(var_names[0]);
 
             if (is_linear) {
                 corrected_profile[element_name] += (reference - measured) * coefficients[1];
             }
             else {
                 corrected_profile[element_name] *= pow(reference / measured, coefficients[1]);
             }
         }
 
         // Apply second correction (typically particle size) if available
         if (var_names.size() > 1 && mlr->Effective(1) && var_names[1] != element_name) {
             double reference = reference_om_size[1];
             double measured = this->at(var_names[1]);
 
             if (is_linear) {
                 corrected_profile[element_name] += (reference - measured) * coefficients[2];
             }
             else {
                 corrected_profile[element_name] *= pow(reference / measured, coefficients[2]);
             }
         }
 
         // Warn about negative corrected values (except for isotopes)
         if (corrected_profile[element_name] < 0 &&
             info.Role != element_information::role::isotope) {
             std::cerr << "Warning: Corrected value for " << element_name
                 << " is negative: " << corrected_profile[element_name] << std::endl;
         }
     }
 
     return corrected_profile;
 }
 
 
 double Elemental_Profile::CalculateDotProduct(const CMBVector& weights) const
 {
     if (size() != weights.getsize()) {
         return -999;
     }
 
     double sum = 0.0;
     int index = 0;
 
     for (const auto& [element_name, concentration] : *this) {
         sum += concentration * weights[index];
         index++;
     }
 
     return sum;
 }
 
 CMBVector Elemental_Profile::SortByConcentration(bool ascending) const
 {
     CMBVector sorted_result;
     vector<string> processed_elements;
 
     for (size_t i = 0; i < size(); i++) {
         double target_value = ascending ? 1e6 : -1e6;
         string target_element;
 
         for (const auto& [element_name, concentration] : *this) {
             // Skip already processed elements
             if (lookup(processed_elements, element_name) != -1) {
                 continue;
             }
 
             // Find min (ascending) or max (descending)
             bool is_better = ascending ? (concentration < target_value)
                 : (concentration > target_value);
 
             if (is_better) {
                 target_element = element_name;
                 target_value = concentration;
             }
         }
 
         processed_elements.push_back(target_element);
         sorted_result.append(target_element, target_value);
     }
 
     return sorted_result;
 }
 
 vector<string> Elemental_Profile::SelectTopElements(int n, bool ascending) const
 {
     CMBVector sorted = SortByConcentration(ascending);
     vector<string> top_elements;
 
     for (int i = 0; i < n && i < sorted.num; i++) {
         top_elements.push_back(sorted.Label(i));
     }
 
     return top_elements;
 }