327 lines
13 KiB
C++
327 lines
13 KiB
C++
//
|
|
// Created by Elie Baier on 23.02.24.
|
|
//
|
|
|
|
#include <cmath>
|
|
#include "iostream"
|
|
#include "string"
|
|
#include "fstream"
|
|
#include "sstream"
|
|
#include "matplot/matplot.h"
|
|
#include <filesystem>
|
|
|
|
using namespace std;
|
|
namespace plt = matplot;
|
|
namespace fs = std::filesystem;
|
|
|
|
int MIN_ZERO_SEPARATION = 100;
|
|
int POWERED_MIN = 3;
|
|
int POWERED_MAX = 3;
|
|
|
|
long double average(vector<long double> a)
|
|
{
|
|
return reduce(a.begin(), a.end(), 0.0) / a.size();
|
|
}
|
|
|
|
vector<long double> findZeroCrossings(const vector<long double>& data, long double absZero) {
|
|
vector<long double> zeroCrossings;
|
|
for (int i = 1; i < data.size(); ++i) {
|
|
if ((data[i - 1] - absZero) * (data[i] - absZero) < 0 || (data[i - 1] - absZero) * (data[i] - absZero) == 0) {
|
|
|
|
// Cleaning the duplicate zero which are +-MIN_ZERO_SEPARATION
|
|
if(zeroCrossings.size() >= 1) {
|
|
if(zeroCrossings[zeroCrossings.size() - 1] - MIN_ZERO_SEPARATION < i && zeroCrossings[zeroCrossings.size() - 1] + MIN_ZERO_SEPARATION > i) {
|
|
// cout << "Found duplicate zero: " << i << endl;
|
|
continue;
|
|
}
|
|
}
|
|
|
|
zeroCrossings.push_back(i);
|
|
}
|
|
}
|
|
return zeroCrossings;
|
|
}
|
|
|
|
pair<long double, bool> poweredMax(const vector<long double>& data, int n) {
|
|
if (data.empty()) {
|
|
return {0.0, false}; // If vector is empty, return false
|
|
}
|
|
|
|
unordered_map<long double, int> freqMap;
|
|
for (long double value : data) {
|
|
freqMap[value]++;
|
|
}
|
|
|
|
long double maxVal = data[0]; // Initialize maxVal with the first element
|
|
int maxFreq = freqMap[maxVal];
|
|
|
|
for (const auto& entry : freqMap) {
|
|
if (entry.second >= n && entry.first > maxVal) {
|
|
maxVal = entry.first;
|
|
maxFreq = entry.second;
|
|
}
|
|
}
|
|
|
|
bool found = maxFreq >= n;
|
|
|
|
return {maxVal, found};
|
|
}
|
|
|
|
std::pair<long double, bool> poweredMin(const std::vector<long double>& data, int n) {
|
|
if (data.empty()) {
|
|
return {0.0, false}; // If vector is empty, return false
|
|
}
|
|
|
|
std::unordered_map<long double, int> freqMap;
|
|
for (long double value : data) {
|
|
freqMap[value]++;
|
|
}
|
|
|
|
long double minVal = std::numeric_limits<long double>::max(); // Initialize minVal with maximum double value
|
|
int minFreq = freqMap[minVal];
|
|
|
|
for (const auto& entry : freqMap) {
|
|
if (entry.second >= n && entry.first < minVal) {
|
|
minVal = entry.first;
|
|
minFreq = entry.second;
|
|
}
|
|
}
|
|
|
|
bool found = minFreq >= n;
|
|
|
|
return {minVal, found};
|
|
}
|
|
|
|
std::vector<long double> slice(const std::vector<long double>& data, int start, int end) {
|
|
if (start < 0 || end >= static_cast<int>(data.size()) || start > end) {
|
|
std::cerr << "Invalid slice indices." << std::endl;
|
|
return {};
|
|
}
|
|
|
|
return {data.begin() + start, data.begin() + end + 1};
|
|
}
|
|
|
|
ifstream openFile(string path) {
|
|
ifstream file;
|
|
file.open(path);
|
|
if(!file.is_open()) {
|
|
std::cerr << "Error opening file." << std::endl;
|
|
}
|
|
return file;
|
|
}
|
|
|
|
int analyseOne(int i, ofstream &out) {
|
|
cout << setprecision(13);
|
|
out << setprecision(13);
|
|
vector<ifstream> channels;
|
|
if(i < 10) {
|
|
channels.push_back(openFile("/Users/eliebaier/Workspace/Physique/TP8/DATA/ALL000" + to_string(i) + "/F000" + to_string(i) + "CH1.CSV"));
|
|
channels.push_back(openFile("/Users/eliebaier/Workspace/Physique/TP8/DATA/ALL000" + to_string(i) + "/F000" + to_string(i) + "CH2.CSV"));
|
|
} else {
|
|
channels.push_back(openFile("/Users/eliebaier/Workspace/Physique/TP8/DATA/ALL00" + to_string(i) + "/F00" + to_string(i) + "CH1.CSV"));
|
|
channels.push_back(openFile("/Users/eliebaier/Workspace/Physique/TP8/DATA/ALL00" + to_string(i) + "/F00" + to_string(i) + "CH2.CSV"));
|
|
}
|
|
|
|
vector<vector<long double>> channelsXData, channelsYData, channelsZeros;
|
|
vector<long double> channelsP2P, channelsAbsoluteZero, channelsTimings, channelsFrequency;
|
|
long double channelOneT = 0;
|
|
|
|
for(int channel = 0; channel < channels.size(); channel++) {
|
|
string line;
|
|
|
|
cout << "========== CHANNEL " << channel << " ==========" << endl;
|
|
cout << "[INFO] Starting data fetching for channel " << channel << endl;
|
|
|
|
// Skipping the TEKTRONICS DEFINITION lines
|
|
cout << "[INFO] Skipping lines used by TEKTRONICS to store device data" << endl;
|
|
for(auto i = 0; i < 18; i++) {
|
|
getline(channels[channel], line);
|
|
}
|
|
|
|
cout << "[INFO] Reading channel data" << endl;
|
|
vector<long double> xs, ys;
|
|
while(getline(channels[channel], line)) {
|
|
// Cleaning the lines
|
|
line.erase(0, 3);
|
|
line.erase(line.size() - 1, line.size());
|
|
|
|
std::string cleanedLine;
|
|
for (char c : line) {
|
|
if (std::isdigit(c) || c == '.' || c == ',') {
|
|
cleanedLine += c;
|
|
}
|
|
}
|
|
|
|
long double x, y;
|
|
istringstream iss(line);
|
|
char comma;
|
|
if (!(iss >> x >> comma >> y)) {
|
|
std::cerr << "Error parsing line: " << line << std::endl;
|
|
continue;
|
|
}
|
|
xs.push_back(x);
|
|
ys.push_back(y);
|
|
}
|
|
channelsXData.push_back(xs);
|
|
channelsYData.push_back(ys);
|
|
cout << "[INFO] Closing channel " << channel << " file" << endl;
|
|
channels[channel].close();
|
|
|
|
// Finding the zeros and finding the average minimum and the average maximum (also finding the reference zero)
|
|
// Finding the absolute max (powered)
|
|
cout << "[INFO] Finding the absolute max and min" << endl;
|
|
double absoluteMax = poweredMax(channelsYData[channel], POWERED_MAX).first;
|
|
cout << "[INFO] Absolute max is " << absoluteMax << endl;
|
|
double absoluteMin = poweredMin(channelsYData[channel], POWERED_MAX).first;
|
|
cout << "[INFO] Absolute min is " << absoluteMin << endl;
|
|
|
|
// Find the reference zero and the offset to the measured one
|
|
channelsAbsoluteZero.push_back((absoluteMax + absoluteMin)/2);
|
|
cout << "[INFO] Reference zero is " << channelsAbsoluteZero[channel] << endl;
|
|
|
|
cout << "[INFO] Finding zeros" << endl;
|
|
channelsZeros.push_back(findZeroCrossings(channelsYData[channel], channelsAbsoluteZero[channel]));
|
|
vector<long double> mins, maxs;
|
|
for(int i = 1; i < channelsZeros[channel].size(); i++) {
|
|
// Checking if it's a min or a max by taking the middle and looking at the sign
|
|
// Add subtraction with channelsAbsoluteZero[channel]
|
|
if(channelsYData[channel][int (channelsZeros[channel][i - 1] + channelsZeros[channel][i])/2] >= 0) {
|
|
// This is a max
|
|
// cout << "Found max between " << xData[zeroCrossings[i - 1]] << " and " << xData[zeroCrossings[i]] << " ";
|
|
// cout << "max is: " << poweredMax(slice(yData, zeroCrossings[i - 1], zeroCrossings[i]), 5).first << endl;
|
|
long double newlyMax = poweredMax(slice(channelsYData[channel], channelsZeros[channel][i - 1], channelsZeros[channel][i]), POWERED_MAX).first;
|
|
cout << "Max found " << newlyMax << endl;
|
|
maxs.push_back(newlyMax);
|
|
} else {
|
|
mins.push_back(poweredMin(slice(channelsYData[channel], channelsZeros[channel][i - 1], channelsZeros[channel][i]), POWERED_MIN).first);
|
|
}
|
|
}
|
|
|
|
vector<long double> frequency;
|
|
for(int j = 0; j + 2 < channelsZeros[channel].size(); j=j+2) {
|
|
// cout << "[INFO] Processing zeros " << j << " and " << j+1 << endl;
|
|
frequency.push_back(((channelsXData[channel][channelsZeros[channel][j + 2]]) - channelsXData[channel][channelsZeros[channel][j]]));
|
|
cout << "[INFO] Distance between zero " << j << " (" << channelsZeros[channel][j] << ") and " << j+2 << " (" << channelsZeros[channel][j + 2] << ") is " << ((channelsXData[channel][channelsZeros[channel][j + 2]]) - channelsXData[channel][channelsZeros[channel][j]]) << " (" << channelsZeros[channel][j+2] - channelsZeros[channel][j] << ")" << endl;
|
|
}
|
|
|
|
long double avgFrequency = average(frequency);
|
|
// avgFrequency *= 2;
|
|
cout << "[INFO] Average peak time for channel " << channel << " is " << avgFrequency << "s" << endl;
|
|
channelsTimings.push_back(avgFrequency);
|
|
cout << "[INFO] Average frequency for channel " << channel << " is " << 1/avgFrequency << "Hz" << endl;
|
|
channelsFrequency.push_back(1/avgFrequency);
|
|
if(channel == 0) {
|
|
channelOneT = avgFrequency;
|
|
// cout << "Set chan 1 freq: " << channelOneT << endl;
|
|
}
|
|
|
|
long double avgMax = 0, avgMin = 0;
|
|
for(auto max : maxs) {
|
|
avgMax += max;
|
|
}
|
|
for(auto min : mins) {
|
|
avgMin += min;
|
|
}
|
|
avgMax /= maxs.size();
|
|
avgMin /= mins.size();
|
|
cout << "[INFO] Average Peak for channel " << channel << " max is " << avgMax << endl;
|
|
cout << "[INFO] Average Peak for channel " << channel << " min is " << avgMin << endl;
|
|
cout << "[INFO] Average Peak-to-Peak for channel " << channel << " is " << abs(avgMin) + abs(avgMax) << endl;
|
|
// Change this to subtraction for cases which are all negative
|
|
channelsP2P.push_back(abs(avgMin) + abs(avgMax));
|
|
}
|
|
|
|
cout << "[INFO] Starting the calculate the phase shift" << endl;
|
|
long double deltas = 0;
|
|
int ndeltas = 0;
|
|
for(int i = 0; i < channelsZeros[0].size(); i++) {
|
|
if(channelsZeros[1].size() <= i) continue;
|
|
long double d = channelsXData[1][channelsZeros[1][i]] - channelsXData[0][channelsZeros[0][i]];
|
|
cout << "[INFO] Phase shift for zero " << i << " is " << ((d) / channelOneT) * (-360) << endl;
|
|
deltas += (d / channelOneT) * (-360);
|
|
ndeltas++;
|
|
}
|
|
deltas /= ndeltas;
|
|
// cout << ndeltas << deltas << channelOneT << endl;
|
|
cout << "[INFO] Average phase shift (CHA0 to CHA1) " << deltas << endl;
|
|
|
|
cout << "[INFO] Writing to CSV" << endl;
|
|
out << i << "," << deltas << "," << channelsP2P[0] << "," << channelsP2P[1] << "," << channelsTimings[0] << "," << channelsTimings[1] << "," << channelsFrequency[0] << "," << channelsFrequency[1] << "\n";
|
|
|
|
plt::tiledlayout(1, 1);
|
|
|
|
auto ax1 = plt::nexttile();
|
|
|
|
plt::plot(channelsXData[0], channelsYData[0]);
|
|
plt::xlabel("Temps [s]");
|
|
plt::ylabel("Tension [V]");
|
|
plt::title("Tension en fonction du temps (Mesure 00" + to_string(i) + ")");
|
|
|
|
plt::hold(matplot::on);
|
|
|
|
plt::plot(channelsXData[1], channelsYData[1])->use_y2(true);
|
|
|
|
plt::grid(matplot::on);
|
|
plt::gca()->minor_grid(true);
|
|
|
|
auto lgd = ::matplot::legend(ax1, {"Générateur", "Bornes de la résistance"});
|
|
lgd->location(matplot::legend::general_alignment::bottomleft);
|
|
lgd->box(false);
|
|
|
|
/* auto ax2 = plt::nexttile();
|
|
plt::plot(ax2, channelsXData[1], channelsYData[1]);
|
|
plt::xlabel(ax2, "Temps [s]");
|
|
plt::ylabel(ax2, "Tension [V]");
|
|
plt::title(ax2, "Tension aux bornes de la résistance en fonction du temps (Mesure 00" + to_string(i) + ")");
|
|
|
|
plt::grid(matplot::on);
|
|
plt::gca()->minor_grid(true); */
|
|
|
|
/* vector<double> zeroCrossings = findZeroCrossings(yData);
|
|
|
|
// Drawing lines on the found zeros
|
|
for(auto zero : zeroCrossings) {
|
|
cout << "Zero: " << zero << endl;
|
|
plt::line(xData[zero], -2, xData[zero], 2);
|
|
} */
|
|
|
|
// Show the plot
|
|
// plt::show();
|
|
|
|
/* cout << "[INFO] Saving svg to " << "/Users/eliebaier/Workspace/Physique/TP8/build/" + to_string(i) + "/ch1.svg" << endl;
|
|
fs::create_directories("./build/" + to_string(i));
|
|
plt::save("/Users/eliebaier/Workspace/Physique/TP8/build/" + to_string(i) + "/ch1.svg"); */
|
|
|
|
/* plt::tiledlayout(2, 1);
|
|
|
|
auto ax2 = plt::nexttile();
|
|
|
|
plt::plot(ax2, channelsXData[1], channelsYData[1]);
|
|
plt::xlabel(ax2, "Temps [s]");
|
|
plt::ylabel(ax2, "Tension [V]");
|
|
plt::title(ax2, "Tension aux bornes de la résistance en fonction du temps (Mesure 00" + to_string(i) + ")"); */
|
|
|
|
/* for(auto zero : channelsZeros[1]) {
|
|
cout << "Zero: " << zero << endl;
|
|
plt::line(channelsXData[1][zero], -2, channelsXData[1][zero], 2);
|
|
} */
|
|
|
|
// plt::show();
|
|
|
|
cout << "[INFO] Saving svg to " << "/Users/eliebaier/Workspace/Physique/TP8/build/" + to_string(i) + "/both-with-label.svg" << endl;
|
|
fs::create_directories("./build/" + to_string(i));
|
|
plt::save("/Users/eliebaier/Workspace/Physique/TP8/build/" + to_string(i) + "/both-with-label.svg");
|
|
|
|
return 0;
|
|
|
|
}
|
|
|
|
int mains() {
|
|
std::ofstream out;
|
|
out.open("/Users/eliebaier/Workspace/Physique/TP8/build/out.csv");
|
|
out << "Mesure, Average Phase Shift, Ch1 Peak-to-Peak, Ch2 Peak-to-Peak, Ch1 Time, Ch2 Time, Ch1 Frequency, Ch2 Frequency" << endl;
|
|
for(int i = 0; i < 53; i++) {
|
|
analyseOne(i, out);
|
|
}
|
|
out.close();
|
|
} |