#include "main.h" #include "mesh.h" #include "filter.h" #include "Settings.h" #include "meshPlotter.h" #include "Plotty.h" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include "FtmKalman.h" #include "mainFtm.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::data.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::data.CurrentPath.NUCs.at(wifi.getAP().getMAC()).ID; float ftm_offset = Settings::data.CurrentPath.NUCs.at(wifi.getAP().getMAC()).ftm_offset; float rssi_pathloss = Settings::data.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::data.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 Stats::Statistics 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]); // ground truth std::vector gtPath = setup.gtPath; Interpolator gtInterpolator = fr.getGroundTruthPath(map, gtPath); Stats::Statistics 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"); // 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); const Point3 srcPath0(9.8, 24.9, 0); // fixed start pos // add shortest-path to destination //const Point3 dst(51, 45, 1.7); //const Point3 dst(25, 45, 0); //NM::NavMeshDijkstra::stamp(mesh, dst); // 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::data.CurrentPath.nucInfo(0).position.xy(), 0.05); plot.addCircle(100000 + 1, Settings::data.CurrentPath.nucInfo(1).position.xy(), 0.05); plot.addCircle(100000 + 2, Settings::data.CurrentPath.nucInfo(2).position.xy(), 0.05); plot.addCircle(100000 + 3, Settings::data.CurrentPath.nucInfo(3).position.xy(), 0.05); 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->kalmanMap = kalmanMap; auto trans = std::make_unique(); //auto trans = std::make_unique(); auto resample = 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; StepDetection sd; PoseDetection pd; TurnDetection td(&pd); RelativePressure relBaro; ActivityDetector act; relBaro.setCalibrationTimeframe( Timestamp::fromMS(5000) ); Timestamp lastTimestamp = Timestamp::fromMS(0); int i = 0; // parse each sensor-value within the offline data for (const Offline::Entry& e : fr.getEntries()) { const Timestamp ts = Timestamp::fromMS(e.ts); if (e.type == Offline::Sensor::WIFI_FTM) { auto ftm = fr.getWifiFtm()[e.idx].data; float ftm_offset = Settings::data.CurrentPath.NUCs.at(ftm.getAP().getMAC()).ftm_offset; float ftmDist = ftm.getFtmDist() + ftm_offset; // in m; plus static offset if (Settings::UseKalman) { auto& kalman = kalmanMap->at(ftm.getAP().getMAC()); float predictDist = kalman.predict(ts, ftmDist); ftm.setFtmDist(predictDist); obs.wifi.insert_or_assign(ftm.getAP().getMAC(), ftm); } else { // MOV AVG if (obs.wifi.count(ftm.getAP().getMAC()) == 0) { obs.wifi.insert_or_assign(ftm.getAP().getMAC(), ftm); } else { auto currFtm = obs.wifi.find(ftm.getAP().getMAC()); float currDist = currFtm->second.getFtmDist(); const float alpha = 0.6; float newDist = alpha * currDist + (1 - alpha) * ftmDist; currFtm->second.setFtmDist(newDist); } } } else if (e.type == Offline::Sensor::WIFI) { //obs.wifi = fr.getWiFiGroupedByTime()[e.idx].data; //ctrl.wifi = fr.getWiFiGroupedByTime()[e.idx].data; } else if (e.type == Offline::Sensor::ACC) { if (sd.add(ts, fr.getAccelerometer()[e.idx].data)) { ++ctrl.numStepsSinceLastEval; } const Offline::TS& _acc = fr.getAccelerometer()[e.idx]; pd.addAccelerometer(ts, _acc.data); //simpleActivity walking / standing act.add(ts, fr.getAccelerometer()[e.idx].data); } else if (e.type == Offline::Sensor::GYRO) { const Offline::TS& _gyr = fr.getGyroscope()[e.idx]; const float delta_gyro = td.addGyroscope(ts, _gyr.data); ctrl.headingChangeSinceLastEval += delta_gyro; } else if (e.type == Offline::Sensor::BARO) { relBaro.add(ts, fr.getBarometer()[e.idx].data); obs.relativePressure = relBaro.getPressureRealtiveToStart(); obs.sigmaPressure = relBaro.getSigma(); //simpleActivity stairs up / down act.add(ts, fr.getBarometer()[e.idx].data); obs.activity = act.get(); } if (ctrl.numStepsSinceLastEval > 0) //if (ts - lastTimestamp >= Timestamp::fromMS(500)) //if (obs.wifi.size() == 4) { obs.currentTime = ts; ctrl.currentTime = ts; // if(ctrl.numStepsSinceLastEval > 0){ // pf.updateTransitionOnly(&ctrl); // } MyState est = pf.update(&ctrl, obs); //pf.updateEvaluationOnly(obs); ctrl.afterEval(); Point3 gtPos = gtInterpolator.get(static_cast(ts.ms())) + Point3(0,0,0.1); lastTimestamp = ts; ctrl.lastEstimate = est.pos.pos; // draw wifi ranges for (auto& ftm : obs.wifi) { int nucid = Settings::data.CurrentPath.NUCs.at(ftm.second.getAP().getMAC()).ID; if (nucid == 1) { Point3 apPos = Settings::data.CurrentPath.NUCs.find(ftm.first)->second.position; //plot.addCircle(nucid, apPos.xy(), ftm.second.getFtmDist()); } } obs.wifi.clear(); //plot //dbg.showParticles(pf.getParticles()); //dbg.setCurPos(est.pos.pos); //dbg.setGT(gtPos); //dbg.addEstimationNode(est.pos.pos); //dbg.addGroundTruthNode(gtPos); //dbg.setTimeInMinute(static_cast(ts.sec()) / 60, static_cast(static_cast(ts.sec())%60)); //dbg.draw(); //plot.printOverview("test"); plot.showParticles(pf.getParticles()); plot.setCurEst(est.pos.pos); plot.setGroundTruth(gtPos); plot.addEstimationNode(est.pos.pos); plot.setActivity((int) act.get()); //plot.splot.getView().setEnabled(false); //plot.splot.getView().setCamera(0, 0); //plot.splot.getView().setEqualXY(true); // plot.plot(); plot.plot(); //plot.closeStream(); std::this_thread::sleep_for(100ms); // error calc // float err_m = gtPos.getDistance(est.pos.pos); // errorStats.add(err_m); // errorFile << ts.ms() << " " << err_m << "\n"; //error calc with penalty for wrong floor double errorFactor = 3.0; Point3 gtPosError = Point3(gtPos.x, gtPos.y, errorFactor * gtPos.z); Point3 estError = Point3(est.pos.pos.x, est.pos.pos.y, errorFactor * est.pos.pos.z); float err_m = gtPosError.getDistance(estError); errorStats.add(err_m); errorFile << ts.ms() << " " << err_m << "\n"; } } // get someting on console std::cout << "Statistical Analysis Filtering: " << std::endl; std::cout << "Median: " << errorStats.getMedian() << " Average: " << errorStats.getAvg() << " Std: " << errorStats.getStdDev() << std::endl; // save the statistical data in file errorFile << "========================================================== \n"; errorFile << "Average of all statistical data: \n"; errorFile << "Median: " << errorStats.getMedian() << "\n"; errorFile << "Average: " << errorStats.getAvg() << "\n"; errorFile << "Standard Deviation: " << errorStats.getStdDev() << "\n"; errorFile << "75 Quantil: " << errorStats.getQuantile(0.75) << "\n"; errorFile.close(); 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); } //mainFtm(argc, argv); //return 0; Stats::Statistics statsAVG; Stats::Statistics statsMedian; Stats::Statistics statsSTD; Stats::Statistics statsQuantil; Stats::Statistics tmp; std::string evaluationName = "prologic/tmp"; 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 < 6; 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(); //for(int i = 0; i < 2; ++i){ // // tmp = run(Settings::data.CurrentPath, 0, evaluationName); // statsMedian.add(tmp.getMedian()); // statsAVG.add(tmp.getAvg()); // statsSTD.add(tmp.getStdDev()); // statsQuantil.add(tmp.getQuantile(0.75)); // std::cout << "Iteration " << i << " completed" << std::endl; //} std::cout << "==========================================================" << std::endl; std::cout << "Average of all statistical data: " << std::endl; std::cout << "Median: " << statsMedian.getAvg() << std::endl; std::cout << "Average: " << statsAVG.getAvg() << std::endl; std::cout << "Standard Deviation: " << statsSTD.getAvg() << std::endl; std::cout << "75 Quantil: " << statsQuantil.getAvg() << std::endl; std::cout << "==========================================================" << std::endl; //EDIT THIS EDIT THIS EDIT THIS EDIT THIS EDIT THIS EDIT THIS EDIT THIS EDIT THIS std::ofstream finalStatisticFile; finalStatisticFile.open (Settings::errorDir + evaluationName + ".csv", std::ios_base::app); finalStatisticFile << "========================================================== \n"; finalStatisticFile << "Average of all statistical data: \n"; finalStatisticFile << "Median: " << statsMedian.getAvg() << "\n"; finalStatisticFile << "Average: " << statsAVG.getAvg() << "\n"; finalStatisticFile << "Standard Deviation: " << statsSTD.getAvg() << "\n"; finalStatisticFile << "75 Quantil: " << statsQuantil.getAvg() << "\n"; finalStatisticFile << "========================================================== \n"; finalStatisticFile.close(); //std::this_thread::sleep_for(std::chrono::seconds(60)); }