#include "main.h" #include "mesh.h" #include "filter.h" #include "Settings.h" #include "meshPlotter.h" #include "Plotty.h" #include "Plotta.h" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include "FtmKalman.h" #include "misc.h" #include using namespace std::chrono_literals; std::vector> getFtmValues(Offline::FileReader& fr, Interpolator& gtInterpolator, const MACAddress nuc) { std::vector> result; for (const Offline::Entry& e : fr.getEntries()) { if (e.type == Offline::Sensor::WIFI_FTM) { const Timestamp ts = Timestamp::fromMS(e.ts); Point3 gtPos = gtInterpolator.get(static_cast(ts.ms())) + Point3(0, 0, 1.3); auto wifi = fr.getWifiFtm()[e.idx].data; if (wifi.getAP().getMAC() == nuc) { Point3 apPos = Settings::CurrentPath.NUCs.find(wifi.getAP().getMAC())->second.position; float apDist = gtPos.getDistance(apPos); float ftmDist = wifi.getFtmDist(); float rssi = wifi.getRSSI(); result.push_back({ apDist, ftmDist, rssi }); } } } return result; } std::pair optimizeFtm(std::vector>& values) { std::vector> error; for (float offset = 0; offset < 10.0f; offset += 0.25) { Stats::Statistics diffs; for (const auto& data : values) { float apDist = std::get<0>(data); float ftmDist = std::get<1>(data); ftmDist += offset; float diff = (apDist - ftmDist); diffs.add(diff); } error.push_back({ offset, diffs.getSquaredSumAvg() }); } auto minElement = std::min_element(error.begin(), error.end(), [](std::pair a, std::pair b) { return a.second < b.second; }); std::cout << "Min ftm offset \t" << minElement->first << "\t" << minElement->second << "\n"; return *minElement; } std::pair optimizeRssi(std::vector>& values) { std::vector> error; for (float pathLoss = 2.0f; pathLoss < 4.0f; pathLoss += 0.125) { Stats::Statistics diffs; for (const auto& data : values) { float apDist = std::get<0>(data); float rssi = std::get<2>(data); float rssiDist = LogDistanceModel::rssiToDistance(-40, pathLoss, rssi); float diff = (apDist - rssiDist); diffs.add(diff); } error.push_back({ pathLoss, diffs.getSquaredSumAvg() }); } auto minElement = std::min_element(error.begin(), error.end(), [](std::pair a, std::pair b) { return a.second < b.second; }); std::cout << "Min path loss \t" << minElement->first << "\t" << minElement->second << "\n"; return *minElement; } void optimizeWifiParameters(Offline::FileReader& fr, Interpolator& gtInterpolator) { int i = 1; for (auto nuc : { Settings::NUC1, Settings::NUC2, Settings::NUC3, Settings::NUC4 }) { auto values = getFtmValues(fr, gtInterpolator, nuc); std::cout << "NUC" << i++ << "\n"; optimizeFtm(values); optimizeRssi(values); } } void exportFtmValues(Offline::FileReader& fr, Interpolator& gtInterpolator) { std::fstream fs; fs.open("test.txt", std::fstream::out); fs << "timestamp;nucid;dist;rssiDist;ftmDist;ftmStdDev;numMeas;numSuccesMeas" << "\n"; for (const Offline::Entry& e : fr.getEntries()) { if (e.type == Offline::Sensor::WIFI_FTM) { const Timestamp ts = Timestamp::fromMS(e.ts); Point3 gtPos = gtInterpolator.get(static_cast(ts.ms())) + Point3(0, 0, 1.3); auto wifi = fr.getWifiFtm()[e.idx].data; int nucid = Settings::CurrentPath.NUCs.at(wifi.getAP().getMAC()).ID; float ftm_offset = Settings::CurrentPath.NUCs.at(wifi.getAP().getMAC()).ftm_offset; float rssi_pathloss = Settings::CurrentPath.NUCs.at(wifi.getAP().getMAC()).rssi_pathloss; float rssiDist = LogDistanceModel::rssiToDistance(-40, rssi_pathloss, wifi.getRSSI()); float ftmDist = wifi.getFtmDist() + ftm_offset; //in m; plus static offset float ftmStdDev = wifi.getFtmDistStd(); int numMeas = wifi.getNumAttemptedMeasurements(); int numSuccessMeas = wifi.getNumSuccessfulMeasurements(); Point3 apPos = Settings::CurrentPath.NUCs.find(wifi.getAP().getMAC())->second.position; float apDist = gtPos.getDistance(apPos); fs << ts.ms() << ";" << nucid << ";" << apDist << ";" << rssiDist << ";" << ftmDist << ";" << ftmStdDev << ";" << numMeas << ";" << numSuccessMeas << "\n"; } } fs.close(); } static float kalman_procNoiseDistStdDev = 1.2f; // standard deviation of distance for process noise static float kalman_procNoiseVelStdDev = 0.1f; // standard deviation of velocity for process noise static CombinedStats run(Settings::DataSetup setup, int walkIdx, std::string folder) { // reading file std::string currDir = std::filesystem::current_path().string(); Floorplan::IndoorMap* map = Floorplan::Reader::readFromFile(setup.map); Offline::FileReader fr(setup.training[walkIdx], setup.HasNanoSecondTimestamps); // ground truth std::vector gtPath = setup.gtPath; Interpolator gtInterpolator = fr.getGroundTruthPath(map, gtPath); CombinedStats errorStats; //calculate distance of path std::vector::InterpolatorEntry> gtEntries = gtInterpolator.getEntries(); double distance = 0; for(int i = 1; i < gtEntries.size(); ++i){ distance += gtEntries[i].value.getDistance(gtEntries[i-1].value); } std::cout << "Distance of Path: " << distance << std::endl; // error file const long int t = static_cast(time(NULL)); auto evalDir = std::filesystem::path(Settings::errorDir); evalDir.append(folder); if (!std::filesystem::exists(evalDir)) { std::filesystem::create_directory(evalDir); } std::ofstream errorFile; errorFile.open (evalDir.string() + "/" + std::to_string(walkIdx) + "_" + std::to_string(t) + ".csv"); // Output dir auto outputDir = std::filesystem::path(Settings::outputDir); outputDir.append(Settings::CurrentPath.name + "_" + std::to_string(walkIdx)); if (!std::filesystem::exists(outputDir)) { std::filesystem::create_directories(outputDir); } // wifi auto kalmanMap = std::make_shared>(); kalmanMap->insert({ Settings::NUC1, Kalman(1, setup.NUCs.at(Settings::NUC1).kalman_measStdDev, kalman_procNoiseDistStdDev, kalman_procNoiseVelStdDev) }); kalmanMap->insert({ Settings::NUC2, Kalman(2, setup.NUCs.at(Settings::NUC2).kalman_measStdDev, kalman_procNoiseDistStdDev, kalman_procNoiseVelStdDev) }); kalmanMap->insert({ Settings::NUC3, Kalman(3, setup.NUCs.at(Settings::NUC3).kalman_measStdDev, kalman_procNoiseDistStdDev, kalman_procNoiseVelStdDev) }); kalmanMap->insert({ Settings::NUC4, Kalman(4, setup.NUCs.at(Settings::NUC4).kalman_measStdDev, kalman_procNoiseDistStdDev, kalman_procNoiseVelStdDev) }); std::cout << "Optimal wifi parameters for " << setup.training[walkIdx] << "\n"; optimizeWifiParameters(fr, gtInterpolator); // mesh NM::NavMeshSettings set; set.maxQuality_m = 0.10; // because of narrow hallways and small rooms reduce min. triangle size (default is 0.2) MyNavMesh mesh; MyNavMeshFactory fac(&mesh, set); fac.build(map); // debug show //MeshPlotter dbg; //dbg.addFloors(map); //dbg.addOutline(map); //dbg.addMesh(mesh); ////dbg.addDijkstra(mesh); //dbg.draw(); Plotty plot(map); plot.buildFloorplan(); plot.setGroundTruth(gtPath); plot.setView(30, 0); // APs Positions plot.addCircle(100000 + 0, Settings::CurrentPath.nucInfo(0).position.xy(), 0.1); plot.addCircle(100000 + 1, Settings::CurrentPath.nucInfo(1).position.xy(), 0.1); plot.addCircle(100000 + 2, Settings::CurrentPath.nucInfo(2).position.xy(), 0.1); plot.addCircle(100000 + 3, Settings::CurrentPath.nucInfo(3).position.xy(), 0.1); plot.plot(); // particle-filter const int numParticles = 5000; //auto init = std::make_unique(&mesh, srcPath0); // known position auto init = std::make_unique(&mesh); // uniform distribution auto eval = std::make_unique(); eval->ftmKalmanFilters = kalmanMap; auto trans = std::make_unique(); //auto trans = std::make_unique(); auto resample = std::make_unique>(); auto estimate = std::make_unique>(); //auto estimate = std::make_unique>(); // setup MyFilter pf(numParticles, std::move(init)); pf.setEvaluation(std::move(eval)); pf.setTransition(std::move(trans)); pf.setResampling(std::move(resample)); pf.setEstimation(std::move(estimate)); pf.setNEffThreshold(0.85); // sensors MyControl ctrl; MyObservation obs; Timestamp lastTimestamp = Timestamp::fromMS(0); std::vector errorValuesFtm, errorValuesRssi; std::vector timestamps; std::vector> gtDistances, ftmDistances, rssiDistances; // distance per AP Plotta::Plotta errorPlot("errorPlot", outputDir.string() + "/errorData.py"); Plotta::Plotta distsPlot("distsPlot", outputDir.string() + "/distances.py"); for (const Offline::Entry& e : fr.getEntries()) { if (e.type != Offline::Sensor::WIFI_FTM) { continue; } // TIME const Timestamp ts = Timestamp::fromMS(e.ts); auto wifiFtm = fr.getWifiFtm()[e.idx].data; obs.ftm.push_back(wifiFtm); if (ts - lastTimestamp >= Timestamp::fromMS(500)) { // Do update step Point2 gtPos = gtInterpolator.get(static_cast(ts.ms())).xy(); plot.setGroundTruth(Point3(gtPos.x, gtPos.y, 0.1)); gtDistances.push_back({ gtPos.getDistance(Settings::CurrentPath.nucInfo(0).position.xy()), gtPos.getDistance(Settings::CurrentPath.nucInfo(1).position.xy()), gtPos.getDistance(Settings::CurrentPath.nucInfo(2).position.xy()), gtPos.getDistance(Settings::CurrentPath.nucInfo(3).position.xy()) }); Point3 estPos; float distErrorFtm = 0; float distErrorRssi = 0; // Run PF obs.currentTime = ts; ctrl.currentTime = ts; MyState est = pf.update(&ctrl, obs); ctrl.afterEval(); lastTimestamp = ts; estPos = est.pos.pos; ctrl.lastEstimate = estPos; plot.setCurEst(Point3(estPos.x, estPos.y, 0.1)); plot.addEstimationNode(Point3(estPos.x, estPos.y, 0.1)); // Error distErrorFtm = gtPos.getDistance(estPos.xy()); errorStats.ftm.add(distErrorFtm); // draw wifi ranges plot.clearDistanceCircles(); for (size_t i = 0; i < obs.ftm.size(); i++) { WiFiMeasurement wifi2 = obs.ftm[i]; Point3 apPos = Settings::CurrentPath.nuc(wifi2.getAP().getMAC()).position; K::GnuplotColor color; switch (Settings::CurrentPath.nuc(wifi2.getAP().getMAC()).ID) { case 1: color = K::GnuplotColor::fromRGB(0, 255, 0); break; case 2: color = K::GnuplotColor::fromRGB(0, 0, 255); break; case 3: color = K::GnuplotColor::fromRGB(255, 255, 0); break; default: color = K::GnuplotColor::fromRGB(255, 0, 0); break; } plot.addDistanceCircle(apPos.xy(), wifi2.getFtmDist(), color); } obs.wifi.clear(); obs.ftm.clear(); errorValuesFtm.push_back(distErrorFtm); errorValuesRssi.push_back(distErrorRssi); timestamps.push_back(ts.ms()); // Error plot errorPlot.add("t", timestamps); errorPlot.add("errorFtm", errorValuesFtm); errorPlot.add("errorRssi", errorValuesRssi); errorPlot.frame(); // Distances plot //distsPlot.add("t", timestamps); //distsPlot.add("gtDists", gtDistances); //distsPlot.add("ftmDists", ftmDistances); //distsPlot.frame(); // Plotting plot.showParticles(pf.getParticles()); plot.setCurEst(estPos); plot.setGroundTruth(Point3(gtPos.x, gtPos.y, 0.1)); plot.addEstimationNode(estPos); //plot.setActivity((int)act.get()); //plot.splot.getView().setEnabled(false); //plot.splot.getView().setCamera(0, 0); //plot.splot.getView().setEqualXY(true); plot.plot(); } } printErrorStats(errorStats); std::ofstream plot_out; plot_out.open(outputDir.string() + "/plot.gp"); plot.clearDistanceCircles(); plot.saveToFile(plot_out); std::ofstream errorStats_out; errorStats_out.open(outputDir.string() + "/error_stats.txt"); printErrorStats(errorStats_out, errorStats); errorPlot.frame(); //system("pause"); return errorStats; } int main(int argc, char** argv) { CmdArguments args(argc, argv); if (args.hasFlag("prob")) { std::cout << "Probabilistic" << "\n"; return mainProp(argc, argv); } else if (args.hasFlag("trilat")) { std::cout << "Trilateration" << "\n"; return mainTrilat(argc, argv); } CombinedStats statsAVG; CombinedStats statsMedian; CombinedStats statsSTD; CombinedStats statsQuantil; CombinedStats tmp; std::string evaluationName = "prologic/tmp"; std::vector setupsToRun = { Settings::data.Path7, Settings::data.Path8, Settings::data.Path9 }; for (Settings::DataSetup setupToRun : setupsToRun) { Settings::CurrentPath = setupToRun; for (size_t walkIdx = 0; walkIdx < 1 /*Settings::CurrentPath.training.size()*/; walkIdx++) { std::cout << "Executing walk " << walkIdx << "\n"; for (int i = 0; i < 1; ++i) { std::cout << "Start of iteration " << i << "\n"; tmp = run(Settings::CurrentPath, walkIdx, evaluationName); statsAVG.ftm.add(tmp.ftm.getAvg()); statsMedian.ftm.add(tmp.ftm.getMedian()); statsSTD.ftm.add(tmp.ftm.getStdDev()); statsQuantil.ftm.add(tmp.ftm.getQuantile(0.75)); statsAVG.rssi.add(tmp.rssi.getAvg()); statsMedian.rssi.add(tmp.rssi.getMedian()); statsSTD.rssi.add(tmp.rssi.getStdDev()); statsQuantil.rssi.add(tmp.rssi.getQuantile(0.75)); std::cout << "Iteration " << i << " completed" << std::endl; } } std::cout << "Results for path " << Settings::CurrentPath.name << std::endl; std::cout << "==========================================================" << std::endl; std::cout << "Average of all statistical data FTM: " << std::endl; std::cout << "Median: " << statsMedian.ftm.getAvg() << std::endl; std::cout << "Average: " << statsAVG.ftm.getAvg() << std::endl; std::cout << "Standard Deviation: " << statsSTD.ftm.getAvg() << std::endl; std::cout << "75 Quantil: " << statsQuantil.ftm.getAvg() << std::endl; std::cout << "==========================================================" << std::endl; std::cout << "==========================================================" << std::endl; std::cout << "Average of all statistical data RSSI: " << std::endl; std::cout << "Median: " << statsMedian.rssi.getAvg() << std::endl; std::cout << "Average: " << statsAVG.rssi.getAvg() << std::endl; std::cout << "Standard Deviation: " << statsSTD.rssi.getAvg() << std::endl; std::cout << "75 Quantil: " << statsQuantil.rssi.getAvg() << std::endl; std::cout << "==========================================================" << std::endl; } //std::vector> error; //std::ofstream error_out; //error_out.open(Settings::errorDir + evaluationName + "_error_path1" + ".csv", std::ios_base::app); //for (kalman_procNoiseDistStdDev = 0.8f; kalman_procNoiseDistStdDev < 1.5f; kalman_procNoiseDistStdDev += 0.1f) //{ // for (kalman_procNoiseVelStdDev = 0.1f; kalman_procNoiseVelStdDev < 0.5f; kalman_procNoiseVelStdDev += 0.1f) // { //for (size_t walkIdx = 0; walkIdx < Settings::data.CurrentPath.training.size(); walkIdx++) //{ // std::cout << "Executing walk " << walkIdx << "\n"; // for (int i = 0; i < 1; ++i) // { // std::cout << "Start of iteration " << i << "\n"; // tmp = run(Settings::data.CurrentPath, walkIdx, evaluationName); // statsMedian.add(tmp.getMedian()); // statsAVG.add(tmp.getAvg()); // statsSTD.add(tmp.getStdDev()); // statsQuantil.add(tmp.getQuantile(0.75)); // std::cout << kalman_procNoiseDistStdDev << " " << kalman_procNoiseVelStdDev << std::endl; // std::cout << "Iteration " << i << " completed" << std::endl; // } //} // error.push_back({{ kalman_procNoiseDistStdDev, kalman_procNoiseVelStdDev, statsAVG.getAvg() }}); // auto minElement = std::min_element(error.begin(), error.end(), [](std::array a, std::array b) { // return a[2] < b[2]; // }); // std::cout << "Current min error " << (*minElement)[2] << "\t Q(0)=\t" << (*minElement)[0] << "\t Q(1)=" << (*minElement)[1] << "\n"; // error_out << kalman_procNoiseDistStdDev << ";" << kalman_procNoiseVelStdDev << ";" << statsAVG.getAvg() << std::endl; // // reset stats // statsAVG.reset(); // statsMedian.reset(); // statsSTD.reset(); // statsQuantil.reset(); // } //} //auto minElement = std::min_element(error.begin(), error.end(), [](std::array a, std::array b) { // return a[2] < b[2]; //}); //std::cout << "Global Min error " << (*minElement)[2] << "\t Q(0)=\t" << (*minElement)[0] << "\t Q(1)=" << (*minElement)[1] << "\n"; //error_out.close(); //std::this_thread::sleep_for(std::chrono::seconds(60)); }