diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
new file mode 100644
index 0000000000000000000000000000000000000000..647b1a72f51067dbda2a154f1cbc131dfec81d6e
--- /dev/null
+++ b/.gitlab-ci.yml
@@ -0,0 +1,9 @@
+# compile jpscore and run light tests
+before_script:
+ - mkdir build && cd build/
+
+build:
+ script:
+ - cmake .. -DBUILD_TESTING=ON
+ - make -j$(nproc)
+ - ctest
\ No newline at end of file
diff --git a/Analysis.cpp b/Analysis.cpp
index a22f30e70fef6b6d0d81a14926cc92c13bf5f871..85542b333285c237bbfdc3675e13610806b06f03 100644
--- a/Analysis.cpp
+++ b/Analysis.cpp
@@ -80,7 +80,7 @@ Analysis::Analysis()
_getProfile = false; // Whether make field analysis or not
_outputGraph = false; // Whether output the data for plot the fundamental diagram each frame
_calcIndividualFD = false; //Adjust whether analyze the individual density and velocity of each pedestrian in stationary state (ALWAYS VORONOI-BASED)
- _vComponent = 'B'; // to mark whether x, y or x and y coordinate are used when calculating the velocity
+ _vComponent = "B"; // to mark whether x, y or x and y coordinate are used when calculating the velocity
_grid_size_X = 0.10; // the size of the grid
_grid_size_Y = 0.10;
@@ -92,6 +92,12 @@ Analysis::Analysis()
_cutRadius=1.0;
_circleEdges=6;
_trajFormat=FileFormat::FORMAT_PLAIN;
+ _isOneDimensional=false;
+ _plotGraph=false;
+ _plotTimeseriesA=false;
+ _plotTimeseriesC=false;
+ _plotTimeseriesD=false;
+
}
Analysis::~Analysis()
@@ -180,11 +186,18 @@ void Analysis::InitArgs(ArgumentParser* args)
_cutByCircle = args->GetIsCutByCircle();
_getProfile = args->GetIsGetProfile();
_outputGraph = args->GetIsOutputGraph();
+ _plotGraph = args->GetIsPlotGraph();
+ _plotTimeseriesA=args->GetIsPlotTimeSeriesA();
+ _plotTimeseriesC=args->GetIsPlotTimeSeriesC();
+ _plotTimeseriesD=args->GetIsPlotTimeSeriesD();
+ _isOneDimensional=args->GetIsOneDimensional();
_calcIndividualFD = args->GetIsIndividualFD();
+ _areaIndividualFD= args->GetAreaIndividualFD();
_vComponent = args->GetVComponent();
_grid_size_X = int(args->GetGridSizeX());
_grid_size_Y = int(args->GetGridSizeY());
_geoPoly = ReadGeometry(args->GetGeometryFilename(), _areaForMethod_D);
+ _geometryFileName=args->GetGeometryFilename();
_projectRootDir=args->GetProjectRootDir();
_trajFormat=args->GetFileFormat();
_cutRadius=args->GetCutRadius();
@@ -211,7 +224,7 @@ std::map<int, polygon_2d> Analysis::ReadGeometry(const std::string& geometryFile
//loop over all areas
for(auto&& area: areas)
{
- //search for the subroom that containst that area
+ //search for the subroom that containst that area
for (auto&& it_room : _building->GetAllRooms())
{
for (auto&& it_sub : it_room.second->GetAllSubRooms())
@@ -281,44 +294,46 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
std::map<int , std::vector<int> > _peds_t=data.GetPedsFrame();
for(int frameNr = 0; frameNr < data.GetNumFrames(); frameNr++ )
{
- vector<int> ids=_peds_t[frameNr];
- vector<int> IdInFrame = data.GetIdInFrame(ids);
- vector<double> XInFrame = data.GetXInFrame(frameNr, ids);
- vector<double> YInFrame = data.GetYInFrame(frameNr, ids);
- for( unsigned int i=0;i<IdInFrame.size();i++)
- {
- bool IsInBuilding=false;
- for (auto&& it_room : _building->GetAllRooms())
- {
- for (auto&& it_sub : it_room.second->GetAllSubRooms())
- {
- SubRoom* subroom = it_sub.second.get();
- if(subroom->IsInSubRoom(Point(XInFrame[i]*CMtoM,YInFrame[i]*CMtoM)))
- {
- IsInBuilding=true;
- break;
- }
- }
- if(IsInBuilding)
- {
- break;
- }
- }
- if(false==IsInBuilding)
- {
- Log->Write("Warning:\tAt %dth frame pedestrian at <x=%.4f, y=%.4f> is not in geometry!", frameNr+data.GetMinFrame(), XInFrame[i]*CMtoM, YInFrame[i]*CMtoM );
- }
- }
+ vector<int> ids=_peds_t[frameNr];
+ vector<int> IdInFrame = data.GetIdInFrame(ids);
+ vector<double> XInFrame = data.GetXInFrame(frameNr, ids);
+ vector<double> YInFrame = data.GetYInFrame(frameNr, ids);
+ for( unsigned int i=0;i<IdInFrame.size();i++)
+ {
+ bool IsInBuilding=false;
+ for (auto&& it_room : _building->GetAllRooms())
+ {
+ for (auto&& it_sub : it_room.second->GetAllSubRooms())
+ {
+ SubRoom* subroom = it_sub.second.get();
+ if(subroom->IsInSubRoom(Point(XInFrame[i]*CMtoM,YInFrame[i]*CMtoM)))
+ {
+ IsInBuilding=true;
+ break;
+ }
+ }
+ if(IsInBuilding)
+ {
+ break;
+ }
+ }
+ if(false==IsInBuilding)
+ {
+ Log->Write("Warning:\tAt %dth frame pedestrian at <x=%.4f, y=%.4f> is not in geometry!", frameNr+data.GetMinFrame(), XInFrame[i]*CMtoM, YInFrame[i]*CMtoM );
+ }
+ }
}
//-----------------------------------------------------------------------------------------------------------------------------------------------------------------
if(_DoesUseMethodA) //Method A
{
- for(unsigned int i=0; i<_areaForMethod_A.size(); i++)
+#pragma omp parallel for
+ for(signed int i=0; i<_areaForMethod_A.size(); i++)
{
Method_A method_A ;
method_A.SetMeasurementArea(_areaForMethod_A[i]);
method_A.SetTimeInterval(_deltaT);
+ method_A.SetPlotTimeSeries(_plotTimeseriesA);
bool result_A=method_A.Process(data,_scriptsLocation);
if(result_A)
{
@@ -333,7 +348,8 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
if(_DoesUseMethodB) //Method_B
{
- for(unsigned int i=0; i<_areaForMethod_B.size(); i++)
+#pragma omp parallel for
+ for(signed int i=0; i<_areaForMethod_B.size(); i++)
{
Method_B method_B;
method_B.SetMeasurementArea(_areaForMethod_B[i]);
@@ -351,7 +367,8 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
if(_DoesUseMethodC) //Method C
{
- for(unsigned int i=0; i<_areaForMethod_C.size(); i++)
+#pragma omp parallel for
+ for(signed int i=0; i<_areaForMethod_C.size(); i++)
{
Method_C method_C;
method_C.SetMeasurementArea(_areaForMethod_C[i]);
@@ -359,6 +376,13 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
if(result_C)
{
Log->Write("INFO:\tSuccess with Method C using measurement area id %d!\n",_areaForMethod_C[i]->_id);
+ if(_plotTimeseriesC)
+ {
+ string parameters_Timeseries="python \""+_scriptsLocation+"/_Plot_timeseries_rho_v.py\" -p \""+ _projectRootDir+VORO_LOCATION + "\" -n "+filename+
+ " -f "+boost::lexical_cast<std::string>(data.GetFps());
+ int res=system(parameters_Timeseries.c_str());
+ Log->Write("INFO:\t time series result: %d ",res);
+ }
}
else
{
@@ -369,14 +393,19 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
if(_DoesUseMethodD) //method_D
{
- for(unsigned int i=0; i<_areaForMethod_D.size(); i++)
+#pragma omp parallel for
+ for(signed int i=0; i<_areaForMethod_D.size(); i++)
{
Method_D method_D;
method_D.SetGeometryPolygon(_geoPoly[_areaForMethod_D[i]->_id]);
+ method_D.SetGeometryFileName(_geometryFileName);
method_D.SetGeometryBoundaries(_lowVertexX, _lowVertexY, _highVertexX, _highVertexY);
method_D.SetGridSize(_grid_size_X, _grid_size_Y);
method_D.SetOutputVoronoiCellData(_outputGraph);
+ method_D.SetPlotVoronoiGraph(_plotGraph);
+ method_D.SetDimensional(_isOneDimensional);
method_D.SetCalculateIndividualFD(_calcIndividualFD);
+ method_D.SetAreaIndividualFD(_areaIndividualFD);
method_D.SetCalculateProfiles(_getProfile);
if(_cutByCircle)
{
@@ -386,7 +415,14 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
bool result_D = method_D.Process(data,_scriptsLocation);
if(result_D)
{
- Log->Write("INFO:\tSuccess with Method D using measurement area id %d!\n",_areaForMethod_D[i]->_id);
+ Log->Write("INFO:\tSuccess with Method D using measurement area id %d!\n",_areaForMethod_D[i]->_id);
+ if(_plotTimeseriesD)
+ {
+ string parameters_Timeseries="python \""+_scriptsLocation+"/_Plot_timeseries_rho_v.py\" -p \""+ _projectRootDir+VORO_LOCATION + "\" -n "+filename+
+ " -f "+boost::lexical_cast<std::string>(data.GetFps());
+ int res=system(parameters_Timeseries.c_str());
+ Log->Write("INFO:\t time series result: %d ",res);
+ }
}
else
{
@@ -394,12 +430,6 @@ int Analysis::RunAnalysis(const string& filename, const string& path)
}
}
}
- if(_DoesUseMethodC || _DoesUseMethodD)
- {
- string parameters_Timeseries="python "+_scriptsLocation+"/_Plot_timeseries_rho_v.py -p \""+ _projectRootDir+VORO_LOCATION + "\" -n "+filename+
- " -f "+boost::lexical_cast<std::string>(data.GetFps());
- system(parameters_Timeseries.c_str());
- }
return 0;
}
diff --git a/Analysis.h b/Analysis.h
index d388ff61c7cc53bfe7c3f62f953d63ab27312f9c..86f4484cea94845a45aa6c0f33d4115557e1b408 100644
--- a/Analysis.h
+++ b/Analysis.h
@@ -129,11 +129,18 @@ private:
double _cutRadius;
int _circleEdges;
bool _getProfile; // Whether make field analysis or not
- bool _outputGraph; // Whether output the data for plot the fundamental diagram each frame
+ bool _outputGraph; // Whether output the data for plot the voronoi diagram each frame
+ bool _plotGraph; // Whether plot the voronoi diagram each frame
+ bool _plotTimeseriesA;
+ bool _plotTimeseriesC;
+ bool _plotTimeseriesD;
+ bool _isOneDimensional;
bool _calcIndividualFD; //Adjust whether analyze the individual density and velocity of each pedestrian in stationary state (ALWAYS VORONOI-BASED)
- char _vComponent; // to mark whether x, y or x and y coordinate are used when calculating the velocity
+ polygon_2d _areaIndividualFD;
+ std::string _vComponent; // to mark whether x, y or x and y coordinate are used when calculating the velocity
std::string _projectRootDir;
std::string _scriptsLocation;
+ std::string _geometryFileName;
FileFormat _trajFormat; // format of the trajectory file
std::vector<MeasurementArea_L*> _areaForMethod_A;
diff --git a/CHANGELOG b/CHANGELOG
new file mode 100644
index 0000000000000000000000000000000000000000..1d0a9c9bfe48338f168c47ef84a5915ade7e40d9
--- /dev/null
+++ b/CHANGELOG
@@ -0,0 +1,101 @@
+# JPSreport v0.9
+
+## Added
+
+## Changed
+
+## Fixed
+
+
+
+# JPSreport v0.8
+
+
+## Added
+
+- A switch is added in the infile for `method_D` to turn off plotting Voronoi diagrams. Now it is possible to only output data for the diagram but not plot figures.
+
+- Switches for plotting time series of density and velocity are added for `method_C` and `method_D` in inifile.
+
+- A switch for plotting N-t diagram is added for `method_A` in inifile.
+
+- An option for analyzing one dimensional trajectory data is added in `method_D`.
+
+- Issue a warning when the voronoi cell cannot be calculated.
+
+- A warning will will be given and the program stops if trajectory for a given pedestrian ID is not continuous.
+
+
+## Changed
+
+- Scripts "_Plot_cell_rho.py" and "_Plot_cell_v.py" are modified. Now the geometry is also plotted when plotting voronoi cells.
+
+- The indicator for velocity component can be specified in trajectory files now (.TXT and .XML)
+
+- Scripts "_Plot_FD.py" is modified!
+
+## Fixed
+
+- Output data file "Folw\_NT\_xxxx.dat" is closed before calling script for plotting N-t diagram.
+
+- A bug relating to transformation of units in `method_B` is fixed.
+
+- Fixed error where all trajectories were co-linear.
+
+- A bug for legend in the script "_Plot_timeseries_rho_v.py" is fixed.
+
+- The case that frame ID and Ped ID in trajectory file are not coutinuous can also be analyzed correctly.
+
+- Now when the given file paths in inifile include blank, it still works on windows system.
+
+- when path of trajectory is not given absolutely, the default location is the same folder with the inifile
+
+
+# JPSreport v0.7
+
+## Added
+
+- Added four demos as examples for using JPSreport
+- Added the option for specifying the location of scripts in configuration file.
+- Embedded python scripts (**\_Plot_N\_t.py**, **\_Plot_timeseries\_rho_v.py**) for plotting N-t diagram (Method A), time series of density/velocity diagram (Method C and D) and Voronoi diagrams (Method D).
+- Added python script (**SteadyState.py**) for automatically detecting steady state of pedestrian flow based on time series of density and velocity. When plotting fundamental diagrams normally only data under steady state are used due to its generality.
+- Added python script (**\_Plot_FD.py**) for plotting fundamenatl diagram based on the detected steady state.
+## Changed
+
+- Changed name of some variables in configuration file:
+
+ **measurementAreas** ---> **measurement_areas**
+
+ **Length_in_movement_direction** ---> **length_in_movement_direction**
+
+ **useXComponent** ---> **use_x_component**
+
+ **useYComponent** ---> **use_y_component**
+
+ **halfFrameNumberToUse** ---> **frame_step**
+
+ **timeInterval** ---> **frame_interval**
+
+ **measurementArea** ---> **measurement_area**
+
+ **outputGraph** ---> **output_graph**
+
+ **individualFDdata** ---> **individual_FD**
+
+ **cutByCircle** ---> **cut_by_circle**
+
+ **getProfile** ---> **profiles**
+
+ **scale_x** ---> **grid_size_x**
+
+ **scale_y** ---> **grid_size_y**
+- Changed the data type of frame rate (fps) from integer to float
+
+- Changed the way for dealing with pedestrian outside geometry. In old version JPSreport stops when some pedestrians are outside geometry but now it continue working by
+removing these pedestrians from the list.
+
+- More than one sub rooms in one geometry can be analysed independently.
+
+## Fixed
+
+- Fixed bug for dealing with obstacles inside geometry.
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 4e8015da664d3b3127f6d9862e415a345bf034d6..551121b1e4280e99ba724424225861b6cb70bc2d 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -71,6 +71,21 @@ if(PROCESSOR_COUNT)
set(CTEST_BUILD_FLAGS "-j${PROCESSOR_COUNT}")
endif(PROCESSOR_COUNT)
+# find the correct OpenMP flag
+FIND_PACKAGE(OpenMP)
+if(OPENMP_FOUND)
+ set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
+ set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} ${OpenMP_CXX_FLAGS}")
+ set(CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} ${OpenMP_EXE_LINKER_FLAGS}")
+else(OPENMP_FOUND)
+ if(CMAKE_CXX_COMPILER_ID STREQUAL "Clang")
+ #set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fopenmp")
+ # somehow find_package(openmp) does not work properly with clang
+ else(CMAKE_CXX_COMPILER_ID STREQUAL "Clang")
+ message( STATUS "Disabling OpenMP support" )
+ endif(CMAKE_CXX_COMPILER_ID STREQUAL "Clang")
+endif(OPENMP_FOUND)
+
# test all cpp-files in Utest
if(BUILD_TESTING)
file(GLOB test_files "${CMAKE_TEST_DIR}/*.cpp")
@@ -135,6 +150,7 @@ set ( header_files
#--------------------
+if(NOT CMAKE_GENERATOR MATCHES "Xcode|Visual Studio")
include(CheckCXXCompilerFlag)
CHECK_CXX_COMPILER_FLAG("-std=c++11" COMPILER_SUPPORTS_CXX11)
CHECK_CXX_COMPILER_FLAG("-std=c++0x" COMPILER_SUPPORTS_CXX0X)
@@ -147,6 +163,7 @@ elseif(COMPILER_SUPPORTS_CXX0X)
else()
message(FATAL_ERROR "The compiler ${CMAKE_CXX_COMPILER} has no c++11 support. Please use a different C++ compiler.")
endif()
+endif()
#---------------------
if("${CMAKE_CXX_COMPILER_ID}" STREQUAL "Clang")
message(STATUS "Using Clang " ${CMAKE_CXX_COMPILER_VERSION})
diff --git a/README.md b/README.md
index 9821d42d6b2258c140b5002b34aaf24fedb6f37a..9c84707e79e6acdd2211183538ee658c7548bd0c 100644
--- a/README.md
+++ b/README.md
@@ -1,58 +1,5 @@
-# JPSreport v0.7
+# JPSreport
+================
-## Running
+Documentation can be found [here](http://jupedsim.github.io/jpsreport/)
-```
-./JPSreport.exe inifile.xml"
-```
-
-From the command line trajectory files mentioned in the inifile will be analysed by using the methods and parameters specified in the `inifile`.
-The output results will be saved in the sub-folder of the folder where the inifile exists.
-
-## Added
-
-- Added four demos as examples for using JPSreport
-- Added the option for specifying the location of scripts in configuration file.
-- Embedded python scripts (**\_Plot_N\_t.py**, **\_Plot_timeseries\_rho_v.py**) for plotting N-t diagram (Method A), time series of density/velocity diagram (Method C and D) and Voronoi diagrams (Method D).
-- Added python script (**SteadyState.py**) for automatically detecting steady state of pedestrian flow based on time series of density and velocity. When plotting fundamental diagrams normally only data under steady state are used due to its generality.
-- Added python script (**\_Plot_FD.py**) for plotting fundamenatl diagram based on the detected steady state.
-
-## Changed
-
-- Changed name of some variables in configuration file:
-
- **measurementAreas** ---> **measurement_areas**
-
- **Length_in_movement_direction** ---> **length_in_movement_direction**
-
- **useXComponent** ---> **use_x_component**
-
- **useYComponent** ---> **use_y_component**
-
- **halfFrameNumberToUse** ---> **frame_step**
-
- **timeInterval** ---> **frame_interval**
-
- **measurementArea** ---> **measurement_area**
-
- **outputGraph** ---> **output_graph**
-
- **individualFDdata** ---> **individual_FD**
-
- **cutByCircle** ---> **cut_by_circle**
-
- **getProfile** ---> **profiles**
-
- **scale_x** ---> **grid_size_x**
-
- **scale_y** ---> **grid_size_y**
-- Changed the data type of frame rate (fps) from integer to float
-
-- Changed the way for dealing with pedestrian outside geometry. In old version JPSreport stops when some pedestrians are outside geometry but now it continue working by
-removing these pedestrians from the list.
-
-- More than one sub rooms in one geometry can be analysed independently.
-
-## Fixed
-
-- Fixed bug for dealing with obstacles inside geometry.
\ No newline at end of file
diff --git a/ReadMe.docx b/ReadMe.docx
deleted file mode 100644
index b3b3bb80990a46c49f71458e5663d237b30bc7e6..0000000000000000000000000000000000000000
Binary files a/ReadMe.docx and /dev/null differ
diff --git a/demos/SteadyState/ReadMe.txt b/demos/SteadyState/ReadMe.txt
index c567866e16f02f68c76a9feb00380d640dabc765..3ba4cfb1219dd64fcdefa6abc6419e987cd738e9 100644
--- a/demos/SteadyState/ReadMe.txt
+++ b/demos/SteadyState/ReadMe.txt
@@ -1 +1,8 @@
-python ../../scripts/SteadyState.py -f ./rho_v_Voronoi_traj_AO_b240.txt_id_1.dat -rs 240 -re 640 -vs 240 -ve 640 -p yes
\ No newline at end of file
+python ../../scripts/SteadyState.py
+-f ./rho_v_Voronoi_traj_AO_b240.txt_id_1.dat # file name with the data
+-a
+-c 0 1 2 # column indexes. 0 for the frame column (first).
+-rs 240 240 # start references
+-re 640 640 # end references
+-xl "t" # xlabel. alternative: "frame"
+-yl "\rho\; / 1/m" "v\; m/s" -p yes # ylabels for every column (optionally with unit)
diff --git a/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png b/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png
index 4824a943296f46e0616c2662ec6bc6c64e156f3d..f8a8121c297b2b683ea3b0cc95ce7efaff8e3e1b 100644
Binary files a/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png and b/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png differ
diff --git a/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt b/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
index e7f0dde0e3d1d2ce0d5030b8f171330357b1d583..674eced85cc88623bd86e8e4edcd5fd515da2ebf 100644
--- a/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
+++ b/demos/SteadyState/SteadyState_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
@@ -1,2 +1,2 @@
-# start end
-121 838
+# start end ratio mean std
+121 838 70.95 3.4702 0.2974
diff --git a/demos/SteadyState/SteadyState_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt b/demos/SteadyState/SteadyState_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
index d96bd39545207aeba2d6de876aa7769cb66679ca..35fe4cae4cf9b999caa9e950ba7ffdd3d7574229 100644
--- a/demos/SteadyState/SteadyState_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
+++ b/demos/SteadyState/SteadyState_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
@@ -1,2 +1,2 @@
-# start end
-202 838
+# start end ratio
+202 838 62.85
diff --git a/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png b/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png
index 6e03c5d6cb3865293e55042ad3bf0718fa67de1c..b1b8b66ff1f42fd07294573128bc0ffcbd8d9393 100644
Binary files a/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png and b/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.png differ
diff --git a/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt b/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
index e80907c104c7753caa725b7d47986457a91e9ab5..c952b47decc3059acefed5c71a9034b53ca9cc08 100644
--- a/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
+++ b/demos/SteadyState/SteadyState_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
@@ -1,2 +1,2 @@
-# start end
-202 938
+# start end ratio mean std
+202 938 72.83 0.7240 0.0323
diff --git a/demos/SteadyState/cusum_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt b/demos/SteadyState/cusum_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
deleted file mode 100644
index 6c9c3252a4d2603ce31a35772ce37eb6c59a5d28..0000000000000000000000000000000000000000
--- a/demos/SteadyState/cusum_rho_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
+++ /dev/null
@@ -1,1014 +0,0 @@
-# frame s
-0 100.0000
-1 100.0000
-2 100.0000
-3 100.0000
-4 100.0000
-5 100.0000
-6 100.0000
-7 100.0000
-8 100.0000
-9 100.0000
-10 100.0000
-11 100.0000
-12 100.0000
-13 100.0000
-14 100.0000
-15 100.0000
-16 100.0000
-17 100.0000
-18 100.0000
-19 100.0000
-20 100.0000
-21 100.0000
-22 100.0000
-23 100.0000
-24 100.0000
-25 100.0000
-26 100.0000
-27 100.0000
-28 100.0000
-29 100.0000
-30 100.0000
-31 100.0000
-32 100.0000
-33 100.0000
-34 100.0000
-35 100.0000
-36 100.0000
-37 100.0000
-38 100.0000
-39 100.0000
-40 100.0000
-41 100.0000
-42 100.0000
-43 100.0000
-44 100.0000
-45 100.0000
-46 100.0000
-47 100.0000
-48 100.0000
-49 100.0000
-50 100.0000
-51 100.0000
-52 100.0000
-53 100.0000
-54 100.0000
-55 100.0000
-56 100.0000
-57 100.0000
-58 100.0000
-59 100.0000
-60 100.0000
-61 100.0000
-62 100.0000
-63 100.0000
-64 100.0000
-65 100.0000
-66 100.0000
-67 99.0000
-68 98.0000
-69 97.0000
-70 96.0000
-71 95.0000
-72 94.0000
-73 93.0000
-74 92.0000
-75 91.0000
-76 90.0000
-77 89.0000
-78 88.0000
-79 87.0000
-80 86.0000
-81 85.0000
-82 84.0000
-83 83.0000
-84 82.0000
-85 81.0000
-86 80.0000
-87 79.0000
-88 78.0000
-89 77.0000
-90 76.0000
-91 75.0000
-92 74.0000
-93 73.0000
-94 72.0000
-95 71.0000
-96 70.0000
-97 69.0000
-98 68.0000
-99 69.0000
-100 70.0000
-101 71.0000
-102 72.0000
-103 73.0000
-104 74.0000
-105 75.0000
-106 76.0000
-107 77.0000
-108 78.0000
-109 79.0000
-110 80.0000
-111 81.0000
-112 82.0000
-113 83.0000
-114 84.0000
-115 85.0000
-116 86.0000
-117 87.0000
-118 88.0000
-119 89.0000
-120 90.0000
-121 91.0000
-122 92.0000
-123 93.0000
-124 94.0000
-125 95.0000
-126 94.0000
-127 93.0000
-128 92.0000
-129 91.0000
-130 90.0000
-131 89.0000
-132 88.0000
-133 87.0000
-134 86.0000
-135 85.0000
-136 84.0000
-137 83.0000
-138 82.0000
-139 81.0000
-140 80.0000
-141 79.0000
-142 78.0000
-143 77.0000
-144 76.0000
-145 75.0000
-146 74.0000
-147 73.0000
-148 72.0000
-149 71.0000
-150 70.0000
-151 69.0000
-152 68.0000
-153 67.0000
-154 66.0000
-155 65.0000
-156 64.0000
-157 63.0000
-158 62.0000
-159 61.0000
-160 60.0000
-161 59.0000
-162 58.0000
-163 57.0000
-164 56.0000
-165 55.0000
-166 54.0000
-167 53.0000
-168 52.0000
-169 51.0000
-170 50.0000
-171 49.0000
-172 48.0000
-173 47.0000
-174 46.0000
-175 45.0000
-176 44.0000
-177 43.0000
-178 42.0000
-179 41.0000
-180 40.0000
-181 39.0000
-182 38.0000
-183 37.0000
-184 36.0000
-185 35.0000
-186 34.0000
-187 33.0000
-188 32.0000
-189 31.0000
-190 30.0000
-191 29.0000
-192 28.0000
-193 27.0000
-194 26.0000
-195 25.0000
-196 24.0000
-197 23.0000
-198 22.0000
-199 21.0000
-200 20.0000
-201 19.0000
-202 18.0000
-203 17.0000
-204 16.0000
-205 15.0000
-206 14.0000
-207 13.0000
-208 12.0000
-209 11.0000
-210 10.0000
-211 9.0000
-212 8.0000
-213 7.0000
-214 6.0000
-215 5.0000
-216 4.0000
-217 3.0000
-218 2.0000
-219 1.0000
-220 0.0000
-221 0.0000
-222 0.0000
-223 0.0000
-224 0.0000
-225 0.0000
-226 0.0000
-227 0.0000
-228 0.0000
-229 0.0000
-230 0.0000
-231 0.0000
-232 0.0000
-233 0.0000
-234 0.0000
-235 0.0000
-236 0.0000
-237 0.0000
-238 0.0000
-239 0.0000
-240 0.0000
-241 0.0000
-242 0.0000
-243 0.0000
-244 0.0000
-245 0.0000
-246 0.0000
-247 0.0000
-248 0.0000
-249 0.0000
-250 0.0000
-251 0.0000
-252 0.0000
-253 0.0000
-254 0.0000
-255 0.0000
-256 0.0000
-257 0.0000
-258 0.0000
-259 0.0000
-260 0.0000
-261 0.0000
-262 0.0000
-263 0.0000
-264 0.0000
-265 0.0000
-266 0.0000
-267 0.0000
-268 0.0000
-269 0.0000
-270 0.0000
-271 0.0000
-272 0.0000
-273 0.0000
-274 0.0000
-275 0.0000
-276 0.0000
-277 0.0000
-278 0.0000
-279 0.0000
-280 0.0000
-281 0.0000
-282 0.0000
-283 0.0000
-284 0.0000
-285 0.0000
-286 0.0000
-287 0.0000
-288 0.0000
-289 0.0000
-290 0.0000
-291 0.0000
-292 0.0000
-293 0.0000
-294 0.0000
-295 0.0000
-296 0.0000
-297 0.0000
-298 0.0000
-299 0.0000
-300 0.0000
-301 0.0000
-302 0.0000
-303 0.0000
-304 0.0000
-305 0.0000
-306 0.0000
-307 0.0000
-308 0.0000
-309 0.0000
-310 0.0000
-311 0.0000
-312 0.0000
-313 0.0000
-314 0.0000
-315 0.0000
-316 0.0000
-317 0.0000
-318 0.0000
-319 0.0000
-320 0.0000
-321 0.0000
-322 0.0000
-323 0.0000
-324 0.0000
-325 0.0000
-326 0.0000
-327 0.0000
-328 0.0000
-329 0.0000
-330 0.0000
-331 0.0000
-332 0.0000
-333 0.0000
-334 0.0000
-335 0.0000
-336 0.0000
-337 0.0000
-338 0.0000
-339 0.0000
-340 0.0000
-341 0.0000
-342 0.0000
-343 0.0000
-344 0.0000
-345 0.0000
-346 0.0000
-347 0.0000
-348 0.0000
-349 0.0000
-350 0.0000
-351 0.0000
-352 0.0000
-353 0.0000
-354 0.0000
-355 0.0000
-356 0.0000
-357 0.0000
-358 0.0000
-359 0.0000
-360 0.0000
-361 0.0000
-362 0.0000
-363 0.0000
-364 0.0000
-365 0.0000
-366 0.0000
-367 0.0000
-368 0.0000
-369 0.0000
-370 0.0000
-371 0.0000
-372 0.0000
-373 0.0000
-374 0.0000
-375 0.0000
-376 0.0000
-377 0.0000
-378 0.0000
-379 0.0000
-380 0.0000
-381 0.0000
-382 0.0000
-383 0.0000
-384 0.0000
-385 0.0000
-386 0.0000
-387 0.0000
-388 0.0000
-389 0.0000
-390 0.0000
-391 0.0000
-392 0.0000
-393 0.0000
-394 0.0000
-395 0.0000
-396 0.0000
-397 0.0000
-398 0.0000
-399 0.0000
-400 0.0000
-401 0.0000
-402 0.0000
-403 0.0000
-404 0.0000
-405 0.0000
-406 0.0000
-407 0.0000
-408 0.0000
-409 0.0000
-410 0.0000
-411 0.0000
-412 0.0000
-413 0.0000
-414 0.0000
-415 0.0000
-416 0.0000
-417 0.0000
-418 0.0000
-419 0.0000
-420 0.0000
-421 0.0000
-422 0.0000
-423 0.0000
-424 0.0000
-425 0.0000
-426 0.0000
-427 0.0000
-428 0.0000
-429 0.0000
-430 0.0000
-431 0.0000
-432 0.0000
-433 0.0000
-434 0.0000
-435 0.0000
-436 0.0000
-437 0.0000
-438 0.0000
-439 0.0000
-440 0.0000
-441 0.0000
-442 0.0000
-443 0.0000
-444 0.0000
-445 0.0000
-446 0.0000
-447 0.0000
-448 0.0000
-449 0.0000
-450 0.0000
-451 0.0000
-452 0.0000
-453 0.0000
-454 0.0000
-455 0.0000
-456 0.0000
-457 0.0000
-458 0.0000
-459 0.0000
-460 0.0000
-461 0.0000
-462 0.0000
-463 0.0000
-464 0.0000
-465 0.0000
-466 0.0000
-467 0.0000
-468 0.0000
-469 0.0000
-470 0.0000
-471 0.0000
-472 0.0000
-473 0.0000
-474 0.0000
-475 0.0000
-476 0.0000
-477 0.0000
-478 0.0000
-479 0.0000
-480 0.0000
-481 0.0000
-482 0.0000
-483 0.0000
-484 0.0000
-485 0.0000
-486 0.0000
-487 0.0000
-488 0.0000
-489 0.0000
-490 0.0000
-491 0.0000
-492 0.0000
-493 0.0000
-494 0.0000
-495 0.0000
-496 0.0000
-497 0.0000
-498 0.0000
-499 0.0000
-500 0.0000
-501 0.0000
-502 0.0000
-503 0.0000
-504 0.0000
-505 0.0000
-506 0.0000
-507 0.0000
-508 0.0000
-509 0.0000
-510 0.0000
-511 0.0000
-512 0.0000
-513 0.0000
-514 0.0000
-515 0.0000
-516 0.0000
-517 0.0000
-518 0.0000
-519 0.0000
-520 0.0000
-521 0.0000
-522 0.0000
-523 0.0000
-524 0.0000
-525 0.0000
-526 0.0000
-527 0.0000
-528 0.0000
-529 0.0000
-530 0.0000
-531 0.0000
-532 0.0000
-533 0.0000
-534 0.0000
-535 0.0000
-536 0.0000
-537 0.0000
-538 0.0000
-539 0.0000
-540 0.0000
-541 0.0000
-542 0.0000
-543 0.0000
-544 0.0000
-545 0.0000
-546 0.0000
-547 0.0000
-548 0.0000
-549 0.0000
-550 0.0000
-551 0.0000
-552 0.0000
-553 0.0000
-554 0.0000
-555 0.0000
-556 0.0000
-557 0.0000
-558 0.0000
-559 0.0000
-560 0.0000
-561 0.0000
-562 0.0000
-563 0.0000
-564 0.0000
-565 0.0000
-566 0.0000
-567 0.0000
-568 0.0000
-569 0.0000
-570 0.0000
-571 0.0000
-572 0.0000
-573 0.0000
-574 0.0000
-575 0.0000
-576 0.0000
-577 0.0000
-578 0.0000
-579 0.0000
-580 0.0000
-581 0.0000
-582 0.0000
-583 0.0000
-584 0.0000
-585 0.0000
-586 0.0000
-587 0.0000
-588 0.0000
-589 0.0000
-590 0.0000
-591 0.0000
-592 0.0000
-593 0.0000
-594 0.0000
-595 0.0000
-596 0.0000
-597 0.0000
-598 0.0000
-599 0.0000
-600 0.0000
-601 0.0000
-602 0.0000
-603 0.0000
-604 0.0000
-605 0.0000
-606 0.0000
-607 0.0000
-608 0.0000
-609 0.0000
-610 0.0000
-611 0.0000
-612 0.0000
-613 0.0000
-614 0.0000
-615 0.0000
-616 0.0000
-617 0.0000
-618 0.0000
-619 0.0000
-620 0.0000
-621 0.0000
-622 0.0000
-623 0.0000
-624 0.0000
-625 0.0000
-626 0.0000
-627 0.0000
-628 0.0000
-629 0.0000
-630 0.0000
-631 0.0000
-632 0.0000
-633 0.0000
-634 0.0000
-635 0.0000
-636 0.0000
-637 0.0000
-638 0.0000
-639 0.0000
-640 0.0000
-641 0.0000
-642 0.0000
-643 0.0000
-644 0.0000
-645 0.0000
-646 0.0000
-647 0.0000
-648 0.0000
-649 0.0000
-650 0.0000
-651 0.0000
-652 0.0000
-653 0.0000
-654 0.0000
-655 0.0000
-656 0.0000
-657 0.0000
-658 0.0000
-659 0.0000
-660 0.0000
-661 0.0000
-662 0.0000
-663 0.0000
-664 0.0000
-665 0.0000
-666 0.0000
-667 0.0000
-668 0.0000
-669 0.0000
-670 0.0000
-671 0.0000
-672 0.0000
-673 0.0000
-674 0.0000
-675 0.0000
-676 0.0000
-677 0.0000
-678 0.0000
-679 0.0000
-680 0.0000
-681 0.0000
-682 0.0000
-683 0.0000
-684 0.0000
-685 0.0000
-686 0.0000
-687 0.0000
-688 0.0000
-689 0.0000
-690 0.0000
-691 0.0000
-692 0.0000
-693 0.0000
-694 0.0000
-695 0.0000
-696 0.0000
-697 0.0000
-698 0.0000
-699 0.0000
-700 0.0000
-701 0.0000
-702 0.0000
-703 0.0000
-704 0.0000
-705 0.0000
-706 0.0000
-707 0.0000
-708 1.0000
-709 2.0000
-710 3.0000
-711 4.0000
-712 5.0000
-713 6.0000
-714 7.0000
-715 8.0000
-716 9.0000
-717 10.0000
-718 11.0000
-719 12.0000
-720 13.0000
-721 14.0000
-722 15.0000
-723 16.0000
-724 17.0000
-725 18.0000
-726 19.0000
-727 20.0000
-728 21.0000
-729 22.0000
-730 23.0000
-731 24.0000
-732 25.0000
-733 26.0000
-734 27.0000
-735 28.0000
-736 29.0000
-737 30.0000
-738 31.0000
-739 32.0000
-740 31.0000
-741 30.0000
-742 29.0000
-743 28.0000
-744 27.0000
-745 26.0000
-746 25.0000
-747 24.0000
-748 23.0000
-749 22.0000
-750 21.0000
-751 20.0000
-752 19.0000
-753 18.0000
-754 17.0000
-755 16.0000
-756 15.0000
-757 14.0000
-758 13.0000
-759 12.0000
-760 11.0000
-761 10.0000
-762 9.0000
-763 8.0000
-764 7.0000
-765 6.0000
-766 5.0000
-767 4.0000
-768 3.0000
-769 2.0000
-770 1.0000
-771 2.0000
-772 3.0000
-773 2.0000
-774 1.0000
-775 0.0000
-776 0.0000
-777 0.0000
-778 0.0000
-779 0.0000
-780 0.0000
-781 0.0000
-782 0.0000
-783 0.0000
-784 0.0000
-785 0.0000
-786 0.0000
-787 0.0000
-788 0.0000
-789 0.0000
-790 0.0000
-791 0.0000
-792 0.0000
-793 0.0000
-794 0.0000
-795 0.0000
-796 0.0000
-797 0.0000
-798 0.0000
-799 0.0000
-800 0.0000
-801 0.0000
-802 0.0000
-803 0.0000
-804 0.0000
-805 0.0000
-806 0.0000
-807 0.0000
-808 0.0000
-809 0.0000
-810 0.0000
-811 0.0000
-812 0.0000
-813 0.0000
-814 0.0000
-815 0.0000
-816 0.0000
-817 0.0000
-818 0.0000
-819 0.0000
-820 0.0000
-821 0.0000
-822 0.0000
-823 0.0000
-824 0.0000
-825 0.0000
-826 0.0000
-827 0.0000
-828 0.0000
-829 0.0000
-830 0.0000
-831 0.0000
-832 0.0000
-833 0.0000
-834 0.0000
-835 0.0000
-836 0.0000
-837 0.0000
-838 0.0000
-839 0.0000
-840 1.0000
-841 2.0000
-842 3.0000
-843 4.0000
-844 5.0000
-845 6.0000
-846 7.0000
-847 8.0000
-848 9.0000
-849 10.0000
-850 11.0000
-851 12.0000
-852 13.0000
-853 14.0000
-854 15.0000
-855 16.0000
-856 17.0000
-857 18.0000
-858 19.0000
-859 20.0000
-860 21.0000
-861 22.0000
-862 23.0000
-863 24.0000
-864 25.0000
-865 26.0000
-866 27.0000
-867 28.0000
-868 29.0000
-869 30.0000
-870 31.0000
-871 32.0000
-872 33.0000
-873 34.0000
-874 35.0000
-875 36.0000
-876 37.0000
-877 38.0000
-878 39.0000
-879 40.0000
-880 41.0000
-881 42.0000
-882 43.0000
-883 44.0000
-884 45.0000
-885 46.0000
-886 47.0000
-887 48.0000
-888 49.0000
-889 50.0000
-890 51.0000
-891 52.0000
-892 53.0000
-893 54.0000
-894 55.0000
-895 56.0000
-896 57.0000
-897 58.0000
-898 59.0000
-899 60.0000
-900 61.0000
-901 62.0000
-902 63.0000
-903 64.0000
-904 65.0000
-905 66.0000
-906 67.0000
-907 68.0000
-908 69.0000
-909 70.0000
-910 71.0000
-911 72.0000
-912 73.0000
-913 74.0000
-914 75.0000
-915 76.0000
-916 77.0000
-917 78.0000
-918 79.0000
-919 80.0000
-920 81.0000
-921 82.0000
-922 83.0000
-923 84.0000
-924 85.0000
-925 86.0000
-926 87.0000
-927 88.0000
-928 89.0000
-929 90.0000
-930 91.0000
-931 92.0000
-932 93.0000
-933 94.0000
-934 95.0000
-935 96.0000
-936 97.0000
-937 98.0000
-938 99.0000
-939 100.0000
-940 100.0000
-941 100.0000
-942 100.0000
-943 100.0000
-944 100.0000
-945 100.0000
-946 100.0000
-947 100.0000
-948 100.0000
-949 100.0000
-950 100.0000
-951 100.0000
-952 100.0000
-953 100.0000
-954 100.0000
-955 100.0000
-956 100.0000
-957 100.0000
-958 100.0000
-959 100.0000
-960 100.0000
-961 100.0000
-962 100.0000
-963 100.0000
-964 100.0000
-965 100.0000
-966 100.0000
-967 100.0000
-968 100.0000
-969 100.0000
-970 100.0000
-971 100.0000
-972 100.0000
-973 100.0000
-974 100.0000
-975 100.0000
-976 100.0000
-977 100.0000
-978 100.0000
-979 100.0000
-980 100.0000
-981 100.0000
-982 100.0000
-983 100.0000
-984 100.0000
-985 100.0000
-986 100.0000
-987 100.0000
-988 100.0000
-989 100.0000
-990 100.0000
-991 100.0000
-992 100.0000
-993 100.0000
-994 100.0000
-995 100.0000
-996 100.0000
-997 100.0000
-998 100.0000
-999 100.0000
-1000 100.0000
-1001 100.0000
-1002 100.0000
-1003 100.0000
-1004 100.0000
-1005 100.0000
-1006 100.0000
-1007 100.0000
-1008 100.0000
-1009 100.0000
-1010 100.0000
-1011 100.0000
-1012 100.0000
diff --git a/demos/SteadyState/cusum_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt b/demos/SteadyState/cusum_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
deleted file mode 100644
index 02d503e053cca5a7f8fad494a58748bf857ce743..0000000000000000000000000000000000000000
--- a/demos/SteadyState/cusum_v_rho_v_Voronoi_traj_AO_b240.txt_id_1.dat.txt
+++ /dev/null
@@ -1,1014 +0,0 @@
-# frame s
-0 100.0000
-1 100.0000
-2 100.0000
-3 100.0000
-4 100.0000
-5 100.0000
-6 100.0000
-7 100.0000
-8 100.0000
-9 100.0000
-10 100.0000
-11 100.0000
-12 100.0000
-13 100.0000
-14 100.0000
-15 100.0000
-16 99.0000
-17 100.0000
-18 100.0000
-19 100.0000
-20 100.0000
-21 100.0000
-22 100.0000
-23 100.0000
-24 100.0000
-25 100.0000
-26 100.0000
-27 100.0000
-28 100.0000
-29 100.0000
-30 100.0000
-31 100.0000
-32 100.0000
-33 100.0000
-34 100.0000
-35 100.0000
-36 100.0000
-37 100.0000
-38 100.0000
-39 100.0000
-40 100.0000
-41 100.0000
-42 100.0000
-43 100.0000
-44 100.0000
-45 100.0000
-46 100.0000
-47 100.0000
-48 100.0000
-49 100.0000
-50 100.0000
-51 100.0000
-52 100.0000
-53 100.0000
-54 100.0000
-55 100.0000
-56 100.0000
-57 100.0000
-58 100.0000
-59 100.0000
-60 100.0000
-61 100.0000
-62 100.0000
-63 100.0000
-64 100.0000
-65 100.0000
-66 100.0000
-67 100.0000
-68 100.0000
-69 100.0000
-70 100.0000
-71 100.0000
-72 100.0000
-73 100.0000
-74 100.0000
-75 100.0000
-76 100.0000
-77 100.0000
-78 100.0000
-79 100.0000
-80 100.0000
-81 100.0000
-82 100.0000
-83 100.0000
-84 100.0000
-85 100.0000
-86 100.0000
-87 100.0000
-88 100.0000
-89 100.0000
-90 100.0000
-91 100.0000
-92 100.0000
-93 100.0000
-94 100.0000
-95 100.0000
-96 100.0000
-97 100.0000
-98 100.0000
-99 100.0000
-100 100.0000
-101 100.0000
-102 100.0000
-103 100.0000
-104 100.0000
-105 100.0000
-106 100.0000
-107 100.0000
-108 100.0000
-109 100.0000
-110 100.0000
-111 100.0000
-112 100.0000
-113 100.0000
-114 100.0000
-115 100.0000
-116 100.0000
-117 100.0000
-118 100.0000
-119 100.0000
-120 100.0000
-121 100.0000
-122 100.0000
-123 100.0000
-124 100.0000
-125 100.0000
-126 100.0000
-127 100.0000
-128 100.0000
-129 100.0000
-130 100.0000
-131 100.0000
-132 100.0000
-133 100.0000
-134 100.0000
-135 100.0000
-136 100.0000
-137 100.0000
-138 100.0000
-139 100.0000
-140 100.0000
-141 100.0000
-142 100.0000
-143 100.0000
-144 100.0000
-145 100.0000
-146 100.0000
-147 100.0000
-148 100.0000
-149 100.0000
-150 100.0000
-151 100.0000
-152 100.0000
-153 100.0000
-154 100.0000
-155 100.0000
-156 100.0000
-157 100.0000
-158 100.0000
-159 100.0000
-160 100.0000
-161 100.0000
-162 100.0000
-163 100.0000
-164 100.0000
-165 100.0000
-166 100.0000
-167 100.0000
-168 100.0000
-169 100.0000
-170 100.0000
-171 100.0000
-172 100.0000
-173 100.0000
-174 100.0000
-175 100.0000
-176 100.0000
-177 100.0000
-178 100.0000
-179 100.0000
-180 100.0000
-181 100.0000
-182 100.0000
-183 100.0000
-184 100.0000
-185 100.0000
-186 100.0000
-187 100.0000
-188 100.0000
-189 100.0000
-190 100.0000
-191 100.0000
-192 100.0000
-193 100.0000
-194 100.0000
-195 100.0000
-196 100.0000
-197 100.0000
-198 100.0000
-199 100.0000
-200 100.0000
-201 100.0000
-202 99.0000
-203 98.0000
-204 97.0000
-205 96.0000
-206 95.0000
-207 94.0000
-208 93.0000
-209 92.0000
-210 91.0000
-211 90.0000
-212 89.0000
-213 88.0000
-214 87.0000
-215 86.0000
-216 85.0000
-217 84.0000
-218 83.0000
-219 82.0000
-220 81.0000
-221 80.0000
-222 79.0000
-223 78.0000
-224 77.0000
-225 76.0000
-226 75.0000
-227 74.0000
-228 73.0000
-229 72.0000
-230 71.0000
-231 70.0000
-232 69.0000
-233 68.0000
-234 67.0000
-235 66.0000
-236 65.0000
-237 64.0000
-238 63.0000
-239 62.0000
-240 61.0000
-241 60.0000
-242 59.0000
-243 58.0000
-244 57.0000
-245 56.0000
-246 55.0000
-247 54.0000
-248 53.0000
-249 52.0000
-250 51.0000
-251 50.0000
-252 49.0000
-253 48.0000
-254 47.0000
-255 46.0000
-256 45.0000
-257 44.0000
-258 43.0000
-259 42.0000
-260 41.0000
-261 40.0000
-262 39.0000
-263 38.0000
-264 37.0000
-265 36.0000
-266 35.0000
-267 34.0000
-268 33.0000
-269 32.0000
-270 31.0000
-271 30.0000
-272 29.0000
-273 28.0000
-274 27.0000
-275 26.0000
-276 25.0000
-277 24.0000
-278 23.0000
-279 22.0000
-280 21.0000
-281 20.0000
-282 19.0000
-283 18.0000
-284 17.0000
-285 16.0000
-286 15.0000
-287 14.0000
-288 13.0000
-289 12.0000
-290 11.0000
-291 10.0000
-292 9.0000
-293 8.0000
-294 7.0000
-295 6.0000
-296 5.0000
-297 4.0000
-298 3.0000
-299 2.0000
-300 1.0000
-301 0.0000
-302 0.0000
-303 0.0000
-304 0.0000
-305 0.0000
-306 0.0000
-307 0.0000
-308 0.0000
-309 0.0000
-310 0.0000
-311 0.0000
-312 0.0000
-313 0.0000
-314 0.0000
-315 0.0000
-316 0.0000
-317 0.0000
-318 0.0000
-319 0.0000
-320 0.0000
-321 0.0000
-322 0.0000
-323 0.0000
-324 0.0000
-325 0.0000
-326 0.0000
-327 0.0000
-328 0.0000
-329 0.0000
-330 0.0000
-331 0.0000
-332 0.0000
-333 0.0000
-334 0.0000
-335 0.0000
-336 0.0000
-337 0.0000
-338 0.0000
-339 0.0000
-340 0.0000
-341 0.0000
-342 0.0000
-343 0.0000
-344 0.0000
-345 0.0000
-346 0.0000
-347 0.0000
-348 0.0000
-349 0.0000
-350 0.0000
-351 0.0000
-352 0.0000
-353 0.0000
-354 0.0000
-355 0.0000
-356 0.0000
-357 0.0000
-358 0.0000
-359 0.0000
-360 0.0000
-361 0.0000
-362 0.0000
-363 0.0000
-364 0.0000
-365 0.0000
-366 0.0000
-367 0.0000
-368 0.0000
-369 0.0000
-370 0.0000
-371 0.0000
-372 0.0000
-373 0.0000
-374 0.0000
-375 0.0000
-376 0.0000
-377 0.0000
-378 0.0000
-379 0.0000
-380 0.0000
-381 0.0000
-382 0.0000
-383 0.0000
-384 0.0000
-385 0.0000
-386 0.0000
-387 0.0000
-388 0.0000
-389 0.0000
-390 0.0000
-391 0.0000
-392 0.0000
-393 0.0000
-394 0.0000
-395 0.0000
-396 0.0000
-397 0.0000
-398 0.0000
-399 0.0000
-400 0.0000
-401 0.0000
-402 0.0000
-403 0.0000
-404 0.0000
-405 0.0000
-406 0.0000
-407 0.0000
-408 0.0000
-409 0.0000
-410 0.0000
-411 0.0000
-412 0.0000
-413 0.0000
-414 0.0000
-415 0.0000
-416 0.0000
-417 0.0000
-418 0.0000
-419 0.0000
-420 0.0000
-421 0.0000
-422 0.0000
-423 0.0000
-424 0.0000
-425 0.0000
-426 0.0000
-427 0.0000
-428 0.0000
-429 0.0000
-430 0.0000
-431 0.0000
-432 0.0000
-433 0.0000
-434 0.0000
-435 0.0000
-436 0.0000
-437 0.0000
-438 0.0000
-439 0.0000
-440 0.0000
-441 0.0000
-442 1.0000
-443 2.0000
-444 1.0000
-445 0.0000
-446 0.0000
-447 0.0000
-448 0.0000
-449 0.0000
-450 0.0000
-451 0.0000
-452 0.0000
-453 0.0000
-454 0.0000
-455 0.0000
-456 0.0000
-457 0.0000
-458 0.0000
-459 0.0000
-460 0.0000
-461 0.0000
-462 0.0000
-463 1.0000
-464 2.0000
-465 3.0000
-466 4.0000
-467 5.0000
-468 6.0000
-469 5.0000
-470 4.0000
-471 3.0000
-472 2.0000
-473 1.0000
-474 0.0000
-475 0.0000
-476 0.0000
-477 0.0000
-478 0.0000
-479 0.0000
-480 0.0000
-481 0.0000
-482 0.0000
-483 0.0000
-484 0.0000
-485 1.0000
-486 2.0000
-487 1.0000
-488 0.0000
-489 0.0000
-490 0.0000
-491 0.0000
-492 0.0000
-493 0.0000
-494 0.0000
-495 0.0000
-496 0.0000
-497 0.0000
-498 0.0000
-499 0.0000
-500 0.0000
-501 0.0000
-502 0.0000
-503 0.0000
-504 0.0000
-505 0.0000
-506 0.0000
-507 0.0000
-508 0.0000
-509 0.0000
-510 0.0000
-511 0.0000
-512 0.0000
-513 0.0000
-514 0.0000
-515 0.0000
-516 0.0000
-517 0.0000
-518 0.0000
-519 0.0000
-520 0.0000
-521 0.0000
-522 0.0000
-523 0.0000
-524 0.0000
-525 0.0000
-526 0.0000
-527 0.0000
-528 0.0000
-529 0.0000
-530 0.0000
-531 0.0000
-532 0.0000
-533 0.0000
-534 0.0000
-535 0.0000
-536 0.0000
-537 0.0000
-538 0.0000
-539 0.0000
-540 0.0000
-541 0.0000
-542 0.0000
-543 0.0000
-544 0.0000
-545 0.0000
-546 0.0000
-547 0.0000
-548 0.0000
-549 0.0000
-550 0.0000
-551 0.0000
-552 0.0000
-553 0.0000
-554 0.0000
-555 0.0000
-556 0.0000
-557 0.0000
-558 0.0000
-559 0.0000
-560 0.0000
-561 0.0000
-562 0.0000
-563 0.0000
-564 0.0000
-565 0.0000
-566 0.0000
-567 0.0000
-568 0.0000
-569 0.0000
-570 0.0000
-571 0.0000
-572 0.0000
-573 0.0000
-574 0.0000
-575 0.0000
-576 0.0000
-577 0.0000
-578 0.0000
-579 0.0000
-580 0.0000
-581 0.0000
-582 0.0000
-583 0.0000
-584 0.0000
-585 0.0000
-586 0.0000
-587 0.0000
-588 0.0000
-589 0.0000
-590 0.0000
-591 0.0000
-592 0.0000
-593 0.0000
-594 0.0000
-595 0.0000
-596 0.0000
-597 0.0000
-598 0.0000
-599 0.0000
-600 0.0000
-601 0.0000
-602 0.0000
-603 0.0000
-604 0.0000
-605 0.0000
-606 0.0000
-607 0.0000
-608 0.0000
-609 0.0000
-610 0.0000
-611 0.0000
-612 0.0000
-613 0.0000
-614 0.0000
-615 0.0000
-616 0.0000
-617 0.0000
-618 0.0000
-619 0.0000
-620 0.0000
-621 0.0000
-622 0.0000
-623 0.0000
-624 0.0000
-625 0.0000
-626 0.0000
-627 0.0000
-628 0.0000
-629 0.0000
-630 0.0000
-631 0.0000
-632 0.0000
-633 0.0000
-634 0.0000
-635 0.0000
-636 0.0000
-637 0.0000
-638 0.0000
-639 0.0000
-640 0.0000
-641 0.0000
-642 0.0000
-643 0.0000
-644 0.0000
-645 0.0000
-646 0.0000
-647 0.0000
-648 0.0000
-649 0.0000
-650 0.0000
-651 0.0000
-652 0.0000
-653 0.0000
-654 0.0000
-655 0.0000
-656 0.0000
-657 0.0000
-658 0.0000
-659 0.0000
-660 0.0000
-661 0.0000
-662 1.0000
-663 2.0000
-664 3.0000
-665 4.0000
-666 5.0000
-667 6.0000
-668 7.0000
-669 8.0000
-670 9.0000
-671 10.0000
-672 11.0000
-673 12.0000
-674 13.0000
-675 14.0000
-676 15.0000
-677 16.0000
-678 17.0000
-679 18.0000
-680 19.0000
-681 20.0000
-682 21.0000
-683 22.0000
-684 23.0000
-685 22.0000
-686 21.0000
-687 20.0000
-688 19.0000
-689 18.0000
-690 17.0000
-691 16.0000
-692 15.0000
-693 14.0000
-694 13.0000
-695 12.0000
-696 11.0000
-697 10.0000
-698 9.0000
-699 8.0000
-700 7.0000
-701 6.0000
-702 5.0000
-703 4.0000
-704 3.0000
-705 2.0000
-706 1.0000
-707 0.0000
-708 0.0000
-709 0.0000
-710 0.0000
-711 0.0000
-712 0.0000
-713 0.0000
-714 0.0000
-715 0.0000
-716 0.0000
-717 0.0000
-718 0.0000
-719 0.0000
-720 0.0000
-721 0.0000
-722 0.0000
-723 0.0000
-724 1.0000
-725 2.0000
-726 1.0000
-727 0.0000
-728 0.0000
-729 0.0000
-730 0.0000
-731 0.0000
-732 0.0000
-733 0.0000
-734 0.0000
-735 0.0000
-736 0.0000
-737 0.0000
-738 0.0000
-739 0.0000
-740 0.0000
-741 0.0000
-742 0.0000
-743 0.0000
-744 0.0000
-745 0.0000
-746 0.0000
-747 0.0000
-748 1.0000
-749 2.0000
-750 3.0000
-751 4.0000
-752 5.0000
-753 6.0000
-754 7.0000
-755 6.0000
-756 5.0000
-757 4.0000
-758 3.0000
-759 2.0000
-760 1.0000
-761 0.0000
-762 0.0000
-763 0.0000
-764 0.0000
-765 0.0000
-766 0.0000
-767 0.0000
-768 0.0000
-769 0.0000
-770 0.0000
-771 0.0000
-772 0.0000
-773 0.0000
-774 0.0000
-775 0.0000
-776 0.0000
-777 0.0000
-778 0.0000
-779 0.0000
-780 0.0000
-781 0.0000
-782 0.0000
-783 0.0000
-784 0.0000
-785 0.0000
-786 0.0000
-787 0.0000
-788 0.0000
-789 0.0000
-790 0.0000
-791 0.0000
-792 0.0000
-793 0.0000
-794 0.0000
-795 0.0000
-796 0.0000
-797 0.0000
-798 0.0000
-799 0.0000
-800 0.0000
-801 0.0000
-802 0.0000
-803 0.0000
-804 0.0000
-805 0.0000
-806 0.0000
-807 0.0000
-808 0.0000
-809 0.0000
-810 0.0000
-811 0.0000
-812 0.0000
-813 0.0000
-814 0.0000
-815 0.0000
-816 0.0000
-817 0.0000
-818 0.0000
-819 0.0000
-820 0.0000
-821 0.0000
-822 0.0000
-823 0.0000
-824 1.0000
-825 2.0000
-826 1.0000
-827 0.0000
-828 0.0000
-829 0.0000
-830 0.0000
-831 0.0000
-832 0.0000
-833 0.0000
-834 0.0000
-835 1.0000
-836 2.0000
-837 3.0000
-838 4.0000
-839 5.0000
-840 6.0000
-841 5.0000
-842 4.0000
-843 3.0000
-844 2.0000
-845 1.0000
-846 0.0000
-847 0.0000
-848 0.0000
-849 0.0000
-850 0.0000
-851 0.0000
-852 0.0000
-853 0.0000
-854 0.0000
-855 0.0000
-856 0.0000
-857 0.0000
-858 0.0000
-859 0.0000
-860 0.0000
-861 1.0000
-862 2.0000
-863 3.0000
-864 4.0000
-865 5.0000
-866 6.0000
-867 7.0000
-868 6.0000
-869 5.0000
-870 4.0000
-871 3.0000
-872 2.0000
-873 1.0000
-874 0.0000
-875 0.0000
-876 0.0000
-877 1.0000
-878 2.0000
-879 3.0000
-880 4.0000
-881 5.0000
-882 6.0000
-883 7.0000
-884 8.0000
-885 9.0000
-886 10.0000
-887 11.0000
-888 12.0000
-889 13.0000
-890 14.0000
-891 15.0000
-892 16.0000
-893 17.0000
-894 18.0000
-895 17.0000
-896 16.0000
-897 15.0000
-898 14.0000
-899 13.0000
-900 12.0000
-901 11.0000
-902 10.0000
-903 9.0000
-904 8.0000
-905 7.0000
-906 6.0000
-907 5.0000
-908 4.0000
-909 3.0000
-910 2.0000
-911 1.0000
-912 0.0000
-913 0.0000
-914 0.0000
-915 0.0000
-916 0.0000
-917 0.0000
-918 0.0000
-919 0.0000
-920 1.0000
-921 2.0000
-922 3.0000
-923 4.0000
-924 5.0000
-925 6.0000
-926 5.0000
-927 4.0000
-928 3.0000
-929 2.0000
-930 1.0000
-931 2.0000
-932 3.0000
-933 4.0000
-934 5.0000
-935 6.0000
-936 7.0000
-937 8.0000
-938 9.0000
-939 10.0000
-940 11.0000
-941 12.0000
-942 13.0000
-943 14.0000
-944 15.0000
-945 16.0000
-946 17.0000
-947 18.0000
-948 19.0000
-949 18.0000
-950 17.0000
-951 16.0000
-952 15.0000
-953 14.0000
-954 13.0000
-955 12.0000
-956 11.0000
-957 10.0000
-958 9.0000
-959 8.0000
-960 7.0000
-961 8.0000
-962 9.0000
-963 10.0000
-964 11.0000
-965 10.0000
-966 9.0000
-967 8.0000
-968 7.0000
-969 6.0000
-970 5.0000
-971 4.0000
-972 3.0000
-973 4.0000
-974 5.0000
-975 6.0000
-976 7.0000
-977 8.0000
-978 9.0000
-979 10.0000
-980 11.0000
-981 12.0000
-982 13.0000
-983 14.0000
-984 15.0000
-985 16.0000
-986 15.0000
-987 14.0000
-988 13.0000
-989 12.0000
-990 11.0000
-991 10.0000
-992 9.0000
-993 8.0000
-994 9.0000
-995 10.0000
-996 11.0000
-997 12.0000
-998 13.0000
-999 14.0000
-1000 15.0000
-1001 16.0000
-1002 17.0000
-1003 18.0000
-1004 19.0000
-1005 20.0000
-1006 21.0000
-1007 22.0000
-1008 23.0000
-1009 24.0000
-1010 25.0000
-1011 26.0000
-1012 27.0000
diff --git a/demos/T-Junction/ini_KO_240_050_240.xml b/demos/T-Junction/ini_KO_240_050_240.xml
index 6dc6e6db715fa75e60dfbf764970e02f9ad0ec72..ae41ff5188abe513a8797da8fab465e79ecdc2cb 100644
--- a/demos/T-Junction/ini_KO_240_050_240.xml
+++ b/demos/T-Junction/ini_KO_240_050_240.xml
@@ -10,6 +10,7 @@
</trajectories>
<!----give relative path based on the location inifile or give the absolute path--->
<scripts location="../../scripts/"/>
+
<measurement_areas unit="m">
<area_B id="1" type="BoundingBox">
@@ -28,9 +29,8 @@
<velocity>
<use_x_component>true</use_x_component>
<use_y_component>true</use_y_component>
-
- <!-- half of the time interval that used to calculate instantaneous velocity
- of ped i [fr] here v_i = (X(t+deltaF) - X(t+deltaF))/(2*deltaF). X is location. -->
+ <!-- The time interval that used to calculate instantaneous velocity
+ of ped i [fr] here v_i = (X(t+frame_step/2) - X(t+frame_step/2))/frame_step. X is location. -->
<frame_step>10</frame_step>
</velocity>
@@ -42,6 +42,7 @@
</frame_interval>
<!-- The coordinate of the line used to calculate the flow and velocity -->
<measurement_area id="2" />
+ <plot_time_series enabled="false"/>
</method_A>
<!-- Method B (Zhang2011a) Vel and Dens based on Tin and Tout -->
@@ -52,13 +53,17 @@
<!-- Method C (Zhang2011a) Classical density and Vel -->
<method_C enabled="true">
<measurement_area id="1" />
+ <plot_time_series enabled="false"/>
</method_C>
<!-- Method D (Zhang2011a) Voronoi density and Vel -->
- <method_D enabled="true" output_graph="true" individual_FD="false">
+ <method_D enabled="true">
<measurement_area id="1" />
- <cut_by_circle enabled="false" radius= "1" edges = "10"/>
- <!-- edges represent the precision of discretization of the circle --->
+ <one_dimensional enabled="false"/>
+ <plot_time_series enabled="false"/>
+ <cut_by_circle enabled="false" radius="1.0" edges="10"/>
+ <output_voronoi_cells enabled="false" plot_graphs="false"/>
+ <individual_FD enabled="false" measurement_area_id="1"/>
<profiles enabled="false" grid_size_x="0.20" grid_size_y="0.20"/>
</method_D>
</JPSreport>
diff --git a/demos/bottleneck/ini_AO_300.xml b/demos/bottleneck/ini_AO_300.xml
index e0151d113b9df3ee4e9478963e2bd2618064b0dd..c351ab1c733fdd84fc815e6e79b43f930024db7c 100644
--- a/demos/bottleneck/ini_AO_300.xml
+++ b/demos/bottleneck/ini_AO_300.xml
@@ -32,8 +32,9 @@
<velocity>
<use_x_component>true</use_x_component>
<use_y_component>true</use_y_component>
- <!-- half of the time interval that used to calculate instantaneous velocity of ped i [fr] here v_i = (X(t+deltaF) - X(t+deltaF))/(2*deltaF). X is location. -->
- <frame_step>5</frame_step>
+ <!-- The time interval that used to calculate instantaneous velocity
+ of ped i [fr] here v_i = (X(t+frame_step/2) - X(t+frame_step/2))/frame_step. X is location. -->
+ <frame_step>10</frame_step>
</velocity>
<!-- Method A (Zhang2011a) Flow and Vel -->
@@ -45,6 +46,7 @@
<!-- The coordinate of the line used to calculate the flow and velocity -->
<measurement_area id="2" />
<measurement_area id="4" />
+ <plot_time_series enabled="true"/>
</method_A>
<!-- Method B (Zhang2011a) Vel and Dens based on Tin and Tout -->
@@ -55,13 +57,17 @@
<!-- Method C (Zhang2011a) Classical density and Vel -->
<method_C enabled="true">
<measurement_area id="1" />
+ <plot_time_series enabled="true"/>
</method_C>
<!-- Method D (Zhang2011a) Voronoi density and Vel -->
- <method_D enabled="true" output_graph="true" individual_FD="true">
- <measurement_area id="1" />
- <cut_by_circle enabled="false" radius= "1" edges = "10"/>
- <!-- edges represent the precision of discretization of the circle --->
+ <method_D enabled="true">
+ <measurement_area id="1" />
+ <one_dimensional enabled="false"/>
+ <plot_time_series enabled="true"/>
+ <cut_by_circle enabled="true" radius="1.0" edges="10"/>
+ <output_voronoi_cells enabled="true" plot_graphs="true"/>
+ <individual_FD enabled="true" measurement_area_id="1"/>
<profiles enabled="false" grid_size_x="0.20" grid_size_y="0.20"/>
</method_D>
</JPSreport>
diff --git a/demos/corridor_JPScore/ini_Aufgabe2.xml b/demos/corridor_JPScore/ini_Aufgabe2.xml
index 00013769e9573c42fbad48be67d375bc11bafe8d..15a25ab2c8f525424d2f8f08532e43e5d468f9b6 100644
--- a/demos/corridor_JPScore/ini_Aufgabe2.xml
+++ b/demos/corridor_JPScore/ini_Aufgabe2.xml
@@ -13,10 +13,10 @@
<measurement_areas unit="m">
<area_B id="1" type="BoundingBox">
- <vertex x="-2.00" y="0.00" />
- <vertex x="-2.00" y="1.50" />
- <vertex x="2.00" y="1.50" />
- <vertex x="2.00" y="0.00" />
+ <vertex x="12.00" y="0.00" />
+ <vertex x="12.00" y="3.50" />
+ <vertex x="15.00" y="3.50" />
+ <vertex x="15.00" y="0.00" />
<length_in_movement_direction distance="4.0" />
</area_B>
<area_L id="2" type="Line">
@@ -28,9 +28,8 @@
<velocity>
<use_x_component>true</use_x_component>
<use_y_component>true</use_y_component>
-
- <!-- half of the time interval that used to calculate instantaneous velocity
- of ped i [fr] here v_i = (X(t+deltaF) - X(t+deltaF))/(2*deltaF). X is location. -->
+ <!-- The time interval that used to calculate instantaneous velocity
+ of ped i [fr] here v_i = (X(t+frame_step/2) - X(t+frame_step/2))/frame_step. X is location. -->
<frame_step>10</frame_step>
</velocity>
@@ -42,6 +41,7 @@
</frame_interval>
<!-- The coordinate of the line used to calculate the flow and velocity -->
<measurement_area id="2" />
+ <plot_time_series enabled="false"/>
</method_A>
<!-- Method B (Zhang2011a) Vel and Dens based on Tin and Tout -->
@@ -50,15 +50,19 @@
</method_B>
<!-- Method C (Zhang2011a) Classical density and Vel -->
- <method_C enabled="false">
+ <method_C enabled="true">
<measurement_area id="1" />
+ <plot_time_series enabled="false"/>
</method_C>
<!-- Method D (Zhang2011a) Voronoi density and Vel -->
- <method_D enabled="false" output_graph="true" individual_FD="false">
+ <method_D enabled="true">
<measurement_area id="1" />
- <cut_by_circle enabled="false" radius= "1" edges = "10"/>
- <!-- edges represent the precision of discretization of the circle --->
+ <one_dimensional enabled="false"/>
+ <plot_time_series enabled="false"/>
+ <cut_by_circle enabled="false" radius="1.0" edges="10"/>
+ <output_voronoi_cells enabled="false" plot_graphs="false"/>
+ <individual_FD enabled="false" measurement_area_id="1"/>
<profiles enabled="false" grid_size_x="0.20" grid_size_y="0.20"/>
</method_D>
</JPSreport>
diff --git a/demos/obstacle/ini_corridor.xml b/demos/obstacle/ini_corridor.xml
index af5dbe510fab371f31701c36982a012afdefcc86..2b174f5a9e972b4991f5c366145885f8965d84f3 100644
--- a/demos/obstacle/ini_corridor.xml
+++ b/demos/obstacle/ini_corridor.xml
@@ -27,8 +27,9 @@
<velocity>
<useXComponent>true</useXComponent>
<useYComponent>false</useYComponent>
- <!-- half of the time interval that used to calculate instantaneous velocity of ped i [fr] here v_i = (X(t+deltaF) - X(t+deltaF))/(2*deltaF). X is location. -->
- <halfFrameNumberToUse>5</halfFrameNumberToUse>
+ <!-- The time interval that used to calculate instantaneous velocity
+ of ped i [fr] here v_i = (X(t+frame_step/2) - X(t+frame_step/2))/frame_step. X is location. -->
+ <halfFrameNumberToUse>10</halfFrameNumberToUse>
</velocity>
<!-- Method A (Zhang2011a) Flow and Vel -->
@@ -39,6 +40,7 @@
</timeInterval>
<!-- The coordinate of the line used to calculate the flow and velocity -->
<measurementArea id="2" />
+ <plot_time_series enabled="false"/>
</method_A>
<!-- Method B (Zhang2011a) Vel and Dens based on Tin and Tout -->
@@ -49,14 +51,18 @@
<!-- Method C (Zhang2011a) Classical density and Vel -->
<method_C enabled="true">
<measurementArea id="1" />
+ <plot_time_series enabled="false"/>
</method_C>
<!-- Method D (Zhang2011a) Voronoi density and Vel -->
- <method_D enabled="true" outputGraph="true" individualFDdata="false">
- <measurementArea id="1" />
- <cutByCircle enabled="false" radius= "1" edges = "10"/>
- <!-- edges represent the precision of discretization of the circle --->
- <getProfile enabled="false" scale_x="0.20" scale_y="0.20"/>
+ <method_D enabled="true">
+ <measurement_area id="1" />
+ <one_dimensional enabled="false"/>
+ <plot_time_series enabled="false"/>
+ <cut_by_circle enabled="false" radius="1.0" edges="10"/>
+ <output_voronoi_cells enabled="true" plot_graphs="false"/>
+ <individual_FD enabled="false" measurement_area_id="1"/>
+ <profiles enabled="false" grid_size_x="0.20" grid_size_y="0.20"/>
</method_D>
</JPSreport>
diff --git a/demos/obstacle/ini_crossing_90_c_6_v5.xml b/demos/obstacle/ini_crossing_90_c_6_v5.xml
index 499e3d4ab972b5ecad502a260bb9b5971ee1313d..0aa60fd4117e86aaa33b4064483d2256194a1823 100644
--- a/demos/obstacle/ini_crossing_90_c_6_v5.xml
+++ b/demos/obstacle/ini_crossing_90_c_6_v5.xml
@@ -28,7 +28,8 @@
<velocity>
<use_x_component>true</use_x_component>
<use_y_component>true</use_y_component>
- <!-- half of the time interval that used to calculate instantaneous velocity of ped i [fr] here v_i = (X(t+deltaF) - X(t+deltaF))/(2*deltaF). X is location. -->
+ <!-- The time interval that used to calculate instantaneous velocity
+ of ped i [fr] here v_i = (X(t+frame_step/2) - X(t+frame_step/2))/frame_step. X is location. -->
<frame_step>5</frame_step>
</velocity>
@@ -40,6 +41,7 @@
</frame_interval>
<!-- The coordinate of the line used to calculate the flow and velocity -->
<measurement_area id="2" />
+ <plot_time_series enabled="false"/>
</method_A>
<!-- Method B (Zhang2011a) Vel and Dens based on Tin and Tout -->
@@ -49,14 +51,18 @@
<!-- Method C (Zhang2011a) Classical density and Vel -->
<method_C enabled="true">
- <measurement_area id="1" />
+ <measurement_area id="1" />
+ <plot_time_series enabled="false"/>
</method_C>
<!-- Method D (Zhang2011a) Voronoi density and Vel -->
- <method_D enabled="true" output_graph="false" individual_FD="false">
+ <method_D enabled="true">
<measurement_area id="1" />
- <cut_by_circle enabled="false" radius= "1" edges = "10"/>
- <!-- edges represent the precision of discretization of the circle --->
+ <one_dimensional enabled="false"/>
+ <plot_time_series enabled="false"/>
+ <cut_by_circle enabled="false" radius="1.0" edges="10"/>
+ <output_voronoi_cells enabled="true" plot_graphs="false"/>
+ <individual_FD enabled="false" measurement_area_id="1"/>
<profiles enabled="false" grid_size_x="0.20" grid_size_y="0.20"/>
</method_D>
</JPSreport>
diff --git a/general/ArgumentParser.cpp b/general/ArgumentParser.cpp
index 173e3c784a56720cb8de9a52276701b7f4fc37f5..f0f5771b709288ea0903289fcfb20fce829494f6 100644
--- a/general/ArgumentParser.cpp
+++ b/general/ArgumentParser.cpp
@@ -27,8 +27,7 @@
**/
-#include <getopt.h>
-#include <unistd.h>
+
#include <cstdio>
#include <cstdlib>
#include <iostream>
@@ -39,6 +38,9 @@
#ifdef _OPENMP
#include <omp.h>
+#else
+#define omp_get_thread_num() 0
+#define omp_get_max_threads() 1
#endif
#include "../tinyxml/tinyxml.h"
@@ -53,47 +55,6 @@ void ArgumentParser::Usage(const std::string file)
Log->Write("Usage: \n");
Log->Write("\t%s input.xml\n",file.c_str());
- /*
- fprintf(stderr,
- "Usage: program options\n\n"
- "With the following options (default values in parenthesis):\n\n"
- " [-t/--trajectory <string>] name of input trajectory file\n"
- " [-I/--input path <filepath>] path of the input file(./Inputfiles/)\n"
- " [-O/--output path <filepath>] path of the input file(./Outputfiles/)\n"
- " [-g/--geometry <string>] path to the geometry file (./Inputfiles/geo.xml)\n"
- " [-m/--measurement area <int>] type of measurement area(1)\n"
- " 1: Bounding box\n"
- " 2: Line\n"
- " [-b/--bounding box <double>] p1.x p1.y p2.x p2.y p3.x p3.y p4.x p4.y (in clockwise)\n"
- " [-d/--moving direction <double>] p1.x p1.y p2.x p2.y \n"
- " [-l/--line <double>] p1.x p1.y p2.x p2.y \n"
-
- fprintf(stderr, "-c --> set cutbycircle=true (false)\n");
- fprintf(stderr, "-a --> set fieldAnalysis=true (false)\n");
- fprintf(stderr, "-g --> set IsOutputGraph=true (false)\n");
- fprintf(stderr, "-v --> set calcIndividualfunddata=true (false)\n");
- fprintf(stderr, "-s scale (3000)\n");
- fprintf(stderr, "-l --> set IsClassicMethod=true (false)\n");
- fprintf(stderr, "-F --> set IsFundamentalTinTout=true (false). density is classical. So IsClassicMethod should be true\n");
- fprintf(stderr, "-V --> set IsFlowVelocity=true (false)\n");
- fprintf(stderr, "-L x1 y1 x2 y2 (0.0, 300.0, 250.0, 300.0)\n");
- fprintf(stderr, "-y beginstationary (700)\n");
- fprintf(stderr, "-Y endstationary (1800)\n");
- fprintf(stderr, "-R Row (65)\n");
- fprintf(stderr, "-C Column (80)\n");
- fprintf(stderr, "-m Meas. Area ax1 (-300)\n");
- fprintf(stderr, "-n Meas. Area ay1 (100)\n");
- fprintf(stderr, "-M Meas. Area ax2 (300)\n");
- fprintf(stderr, "-N Meas. Area ay2 (200)\n");
- fprintf(stderr, "-o Outputfile (result.dat)\n");
- fprintf(stderr, "-O goes in the name of the polygons, speed and point files (dummy)\n");
- fprintf(stderr, "-d --> set use_Vxy false (true)\n");
- fprintf(stderr, "-e --> set use_Vy true (false)\n");
- fprintf(stderr, "-k --> set use_Vx true (false)\n");
- fprintf(stderr, "-p fps (10)\n");
- fprintf(stderr, "-G cor_x cor_y length width (corridor)\n");
-
- */
exit(EXIT_SUCCESS);
}
@@ -102,7 +63,7 @@ ArgumentParser::ArgumentParser()
// Default parameter values
_geometryFileName = "geo.xml";
- _vComponent = 'B';
+ _vComponent = "B";
_isMethodA = false;
_timeIntervalA = 160;
_delatTVInst = 5;
@@ -111,6 +72,11 @@ ArgumentParser::ArgumentParser()
_isMethodD = false;
_isCutByCircle = false;
_isOutputGraph= false;
+ _isPlotGraph= false;
+ _isPlotTimeSeriesA=false;
+ _isPlotTimeSeriesD=false;
+ _isPlotTimeSeriesC=false;
+ _isOneDimensional=false;
_isIndividualFD = false;
_isGetProfile =false;
_steadyStart =100;
@@ -181,7 +147,7 @@ const vector<string>& ArgumentParser::GetTrajectoriesFiles() const
return _trajectoriesFiles;
}
-const string& ArgumentParser::ArgumentParser::GetProjectRootDir() const
+const string& ArgumentParser::GetProjectRootDir() const
{
return _projectRootDir;
}
@@ -219,11 +185,11 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
}
//geometry
- if(xMainNode->FirstChild("geometry"))
- {
- _geometryFileName=_projectRootDir+xMainNode->FirstChildElement("geometry")->Attribute("file");
- Log->Write("INFO: \tGeometry File is: <"+_geometryFileName+">");
- }
+ if(xMainNode->FirstChild("geometry"))
+ {
+ _geometryFileName=_projectRootDir+xMainNode->FirstChildElement("geometry")->Attribute("file");
+ Log->Write("INFO: \tGeometry File is: <"+_geometryFileName+">");
+ }
//trajectories
TiXmlNode* xTrajectories = xMainNode->FirstChild("trajectories");
@@ -275,82 +241,103 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
if (xTrajectories->FirstChildElement("path"))
{
- _trajectoriesLocation = xTrajectories->FirstChildElement("path")->Attribute("location");
- if(_trajectoriesLocation.empty())
- _trajectoriesLocation="./";
-
+ if(xTrajectories->FirstChildElement("path")->Attribute("location"))
+ {
+ _trajectoriesLocation = xTrajectories->FirstChildElement("path")->Attribute("location");
+ }
//hack to find if it is an absolute path
// ignore the project root in this case
if ( (boost::algorithm::contains(_trajectoriesLocation,":")==false) && //windows
(boost::algorithm::starts_with(_trajectoriesLocation,"/") ==false)) //linux
// &&() osx
{
- _trajectoriesLocation=_projectRootDir+_trajectoriesLocation;
+ Log->Write(_trajectoriesLocation);
+ _trajectoriesLocation=_projectRootDir+_trajectoriesLocation;
}
+ }
+ else
+ {
+ _trajectoriesLocation=_projectRootDir;
+ }
- // in the case no file was specified, collect all xml files in the specified directory
- if(_trajectoriesFiles.empty())
+ Log->Write("INFO: \tInput directory for loading trajectory is:\t<"+ (_trajectoriesLocation)+">");
+
+ // in the case no file was specified, collect all files in the specified directory
+ if(_trajectoriesFiles.empty())
+ {
+ DIR *dir;
+ struct dirent *ent;
+ if ((dir = opendir (_trajectoriesLocation.c_str())) != NULL)
{
- DIR *dir;
- struct dirent *ent;
- if ((dir = opendir (_trajectoriesLocation.c_str())) != NULL)
+ /* print all the files and directories within directory */
+ while ((ent = readdir (dir)) != NULL)
{
- /* print all the files and directories within directory */
- while ((ent = readdir (dir)) != NULL)
+ string filename=ent->d_name;
+
+ if (boost::algorithm::ends_with(filename, fmt))
+ //if (filename.find(fmt)!=std::string::npos)
{
- string filename=ent->d_name;
-
- if (boost::algorithm::ends_with(filename, fmt))
- //if (filename.find(fmt)!=std::string::npos)
- {
- //_trajectoriesFiles.push_back(_projectRootDir+filename);
- _trajectoriesFiles.push_back(filename);
- Log->Write("INFO: \tInput trajectory file is\t<"+ (filename)+">");
- }
+ //_trajectoriesFiles.push_back(_projectRootDir+filename);
+ _trajectoriesFiles.push_back(filename);
+ Log->Write("INFO: \tInput trajectory file is\t<"+ (filename)+">");
}
- closedir (dir);
- }
- else
- {
- /* could not open directory */
- Log->Write("ERROR: \tcould not open the directory <"+_trajectoriesLocation+">");
- return false;
}
+ closedir (dir);
+ }
+ else
+ {
+ /* could not open directory */
+ Log->Write("ERROR: \tcould not open the directory <"+_trajectoriesLocation+">");
+ return false;
}
}
- Log->Write("INFO: \tInput directory for loading trajectory is:\t<"+ (_trajectoriesLocation)+">");
+
+ }
+
+ //max CPU
+ if(xMainNode->FirstChild("num_threads"))
+ {
+ TiXmlNode* numthreads = xMainNode->FirstChild("num_threads")->FirstChild();
+ if(numthreads)
+ {
+#ifdef _OPENMP
+ omp_set_num_threads(xmltoi(numthreads->Value(),omp_get_max_threads()));
+#endif
+
+ }
+ Log->Write("INFO: \t Using <%d> threads", omp_get_max_threads());
}
//scripts
- if(xMainNode->FirstChild("scripts"))
- {
- _scriptsLocation=xMainNode->FirstChildElement("scripts")->Attribute("location");
- if(_scriptsLocation.empty())
- {
- _scriptsLocation="./";
- }
- if ( (boost::algorithm::contains(_scriptsLocation,":")==false) && //windows
- (boost::algorithm::starts_with(_scriptsLocation,"/") ==false)) //linux
- // &&() osx
- {
- _scriptsLocation=_projectRootDir+_scriptsLocation;
- }
- if (opendir (_scriptsLocation.c_str()) == NULL)
- {
- /* could not open directory */
- Log->Write("ERROR: \tcould not open the directory <"+_scriptsLocation+">");
- return false;
- }
- else
- {
- Log->Write("INFO: \tInput directory for loading scripts is:\t<"+_scriptsLocation+">");
- }
- }
+ if(xMainNode->FirstChild("scripts"))
+ {
+ _scriptsLocation=xMainNode->FirstChildElement("scripts")->Attribute("location");
+ if(_scriptsLocation.empty())
+ {
+ _scriptsLocation="./";
+ }
+ if ( (boost::algorithm::contains(_scriptsLocation,":")==false) && //windows
+ (boost::algorithm::starts_with(_scriptsLocation,"/") ==false)) //linux
+ // &&() osx
+ {
+ _scriptsLocation=_projectRootDir+_scriptsLocation;
+ }
+ if (opendir (_scriptsLocation.c_str()) == NULL)
+ {
+ /* could not open directory */
+ Log->Write("ERROR: \tcould not open the directory <"+_scriptsLocation+">");
+ return false;
+ }
+ else
+ {
+ Log->Write("INFO: \tInput directory for loading scripts is:\t<"+_scriptsLocation+">");
+ }
+ }
//measurement area
- if(xMainNode->FirstChild("measurement_areas")) {
-
+ if(xMainNode->FirstChild("measurement_areas"))
+ {
string unit = xMainNode->FirstChildElement("measurement_areas")->Attribute("unit");
if(unit!="m")
{
@@ -405,31 +392,70 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
//instantaneous velocity
TiXmlNode* xVelocity=xMainNode->FirstChild("velocity");
- if(xVelocity) {
- string UseXComponent = xVelocity->FirstChildElement("use_x_component")->GetText();
- string UseYComponent = xVelocity->FirstChildElement("use_y_component")->GetText();
- string FrameSteps = xVelocity->FirstChildElement("frame_step")->GetText();
-
- _delatTVInst = atof(FrameSteps.c_str())/2.0;
- if(UseXComponent == "true"&&UseYComponent == "false")
+ if(xVelocity)
+ {
+ string FrameSteps = xVelocity->FirstChildElement("frame_step")->GetText();
+ _delatTVInst = atof(FrameSteps.c_str())/2.0;
+ TiXmlNode* xVx=xVelocity->FirstChildElement("use_x_component");
+ TiXmlNode* xVy=xVelocity->FirstChildElement("use_y_component");
+ //decide which component used in velocity calculation
+ if(xVx && xVy)
{
- _vComponent = 'X';
- Log->Write("INFO: \tOnly x-component coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ string UseXComponent = xVelocity->FirstChildElement("use_x_component")->GetText();
+ string UseYComponent = xVelocity->FirstChildElement("use_y_component")->GetText();
+ if(UseXComponent == "true"&&UseYComponent == "false")
+ {
+ _vComponent = "X";
+ Log->Write("INFO: \tOnly x-component coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
+ else if(UseXComponent == "false"&&UseYComponent == "true")
+ {
+ _vComponent = "Y";
+ Log->Write("INFO: \tOnly y-component coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
+ else if(UseXComponent == "true"&&UseYComponent == "true")
+ {
+ _vComponent = "B"; // both components
+ Log->Write("INFO: \tBoth x and y-component of coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
+ else if(UseXComponent == "false"&&UseYComponent == "false")
+ {
+ _vComponent = "F";
+ Log->Write("INFO: \tThe component defined in the trajectory file will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
}
- else if(UseXComponent == "false"&&UseYComponent == "true")
+ else if(xVx && !xVy)
{
- _vComponent = 'Y';
- Log->Write("INFO: \tOnly y-component coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ string UseXComponent = xVelocity->FirstChildElement("use_x_component")->GetText();
+ if(UseXComponent == "true")
+ {
+ _vComponent = "X";
+ Log->Write("INFO: \tOnly x-component coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
+ else if(UseXComponent == "false")
+ {
+ _vComponent = "F";
+ Log->Write("INFO: \tThe component defined in the trajectory file will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
}
- else if(UseXComponent == "true"&&UseYComponent == "true")
+ else if(!xVx && xVy)
{
- _vComponent = 'B'; // both components
- Log->Write("INFO: \tBoth x and y-component of coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ string UseYComponent = xVelocity->FirstChildElement("use_y_component")->GetText();
+ if(UseYComponent == "true")
+ {
+ _vComponent = "Y";
+ Log->Write("INFO: \tOnly y-component coordinates will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
+ else if(UseYComponent == "false")
+ {
+ _vComponent = "F";
+ Log->Write("INFO: \tThe component defined in the trajectory file will be used to calculate instantaneous velocity over <"+FrameSteps+" frames>" );
+ }
}
else
{
- Log->Write("Error: \tType of velocity is not selected, please check it !!! " );
- return false;
+ _vComponent = "F";
+ Log->Write("INFO: \tThe component defined in the trajectory file will be used to calculate instantaneous velocity over <" + FrameSteps + " frames>");
}
}
@@ -449,6 +475,14 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
_areaIDforMethodA.push_back(xmltoi(xMeasurementArea->Attribute("id")));
Log->Write("INFO: \tMeasurement area id <%d> will be used for analysis", xmltoi(xMeasurementArea->Attribute("id")));
}
+ if(xMethod_A->FirstChildElement("plot_time_series"))
+ {
+ if ( string(xMethod_A->FirstChildElement("plot_time_series")->Attribute("enabled"))=="true")
+ {
+ _isPlotTimeSeriesA=true;
+ Log->Write("INFO: \tThe Time series N-t measured with Method A will be plotted!!");
+ }
+ }
}
}
// method B
@@ -479,6 +513,14 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
_areaIDforMethodC.push_back(xmltoi(xMeasurementArea->Attribute("id")));
Log->Write("INFO: \tMeasurement area id <%d> will be used for analysis", xmltoi(xMeasurementArea->Attribute("id")));
}
+ if (xMethod_C->FirstChildElement("plot_time_series"))
+ {
+ if ( string(xMethod_C->FirstChildElement("plot_time_series")->Attribute("enabled"))=="true")
+ {
+ _isPlotTimeSeriesC=true;
+ Log->Write("INFO: \tThe Time series measured with Method C will be plotted!!");
+ }
+ }
}
// method D
TiXmlElement* xMethod_D=xMainNode->FirstChildElement("method_D");
@@ -487,22 +529,30 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
{
_isMethodD = true;
Log->Write("INFO: \tMethod D is selected" );
- _isOutputGraph = (string(xMethod_D->Attribute("output_graph")) == "true");
- if(_isOutputGraph)
- {
- Log->Write("INFO: \tVoronoi graph is asked to output" );
- }
- _isIndividualFD = (string(xMethod_D->Attribute("individual_FD")) == "true");
- if(_isIndividualFD)
- {
- Log->Write("INFO: \tIndividual fundamental diagram data will be calculated" );
- }
+
for(TiXmlElement* xMeasurementArea=xMainNode->FirstChildElement("method_D")->FirstChildElement("measurement_area");
xMeasurementArea; xMeasurementArea = xMeasurementArea->NextSiblingElement("measurement_area"))
{
_areaIDforMethodD.push_back(xmltoi(xMeasurementArea->Attribute("id")));
Log->Write("INFO: \tMeasurement area id <%d> will be used for analysis", xmltoi(xMeasurementArea->Attribute("id")));
}
+ if (xMethod_D->FirstChildElement("one_dimensional"))
+ {
+ if ( string(xMethod_D->FirstChildElement("one_dimensional")->Attribute("enabled"))=="true")
+ {
+ _isOneDimensional=true;
+ Log->Write("INFO: \tThe data will be analyzed with one dimensional way!!");
+ }
+ }
+
+ if (xMethod_D->FirstChildElement("plot_time_series"))
+ {
+ if ( string(xMethod_D->FirstChildElement("plot_time_series")->Attribute("enabled"))=="true")
+ {
+ _isPlotTimeSeriesD=true;
+ Log->Write("INFO: \tThe Time series measured with Method D will be plotted!!");
+ }
+ }
if ( xMethod_D->FirstChildElement("cut_by_circle"))
{
if ( string(xMethod_D->FirstChildElement("cut_by_circle")->Attribute("enabled"))=="true")
@@ -515,6 +565,34 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
}
}
+ if ( xMethod_D->FirstChildElement("output_voronoi_cells"))
+ {
+ if ( string(xMethod_D->FirstChildElement("output_voronoi_cells")->Attribute("enabled"))=="true")
+ {
+ _isOutputGraph=true;
+ Log->Write("INFO: \tData of voronoi diagram is asked to output" );
+ if(string(xMethod_D->FirstChildElement("output_voronoi_cells")->Attribute("plot_graphs"))=="true")
+ {
+ _isPlotGraph=true;
+ if(_isPlotGraph)
+ {
+ Log->Write("INFO: \tGraph of voronoi diagram will be plotted" );
+ }
+ }
+ }
+ }
+
+ if ( xMethod_D->FirstChildElement("individual_FD"))
+ {
+ if ( string(xMethod_D->FirstChildElement("individual_FD")->Attribute("enabled"))=="true")
+ {
+ _isIndividualFD=true;
+ int areaId=xmltoi(xMethod_D->FirstChildElement("individual_FD")->Attribute("measurement_area_id"));
+ _areaIndividualFD=((MeasurementArea_B *)(_measurementAreas[areaId]))->_poly;
+ Log->Write("INFO: \tIndividual fundamental diagram data will be calculated over the measurement area with id <%d>! ", areaId);
+ }
+ }
+
if ( xMethod_D->FirstChildElement("steadyState"))
{
_steadyStart =xmltof(xMethod_D->FirstChildElement("steadyState")->Attribute("start"));
@@ -523,6 +601,7 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
}
if(xMethod_D->FirstChildElement("profiles"))
+ {
if ( string(xMethod_D->FirstChildElement("profiles")->Attribute("enabled"))=="true")
{
_isGetProfile = true;
@@ -531,6 +610,7 @@ bool ArgumentParser::ParseIniFile(const string& inifile)
Log->Write("INFO: \tProfiles will be calculated" );
Log->Write("INFO: \tThe discretized grid size in x, y direction is: < %f >m by < %f >m ",_grid_size_X*CMtoM, _grid_size_Y*CMtoM);
}
+ }
}
}
Log->Write("INFO: \tFinish parsing inifile");
@@ -573,7 +653,7 @@ const string& ArgumentParser::GetTrajectoriesFilename() const
return _trajectoriesFilename;
}
-char ArgumentParser::GetVComponent() const
+std::string ArgumentParser::GetVComponent() const
{
return _vComponent;
}
@@ -628,11 +708,41 @@ bool ArgumentParser::GetIsOutputGraph() const
return _isOutputGraph;
}
+bool ArgumentParser::GetIsPlotGraph() const
+{
+ return _isPlotGraph;
+}
+
+bool ArgumentParser::GetIsPlotTimeSeriesA() const
+{
+ return _isPlotTimeSeriesA;
+}
+
+bool ArgumentParser::GetIsPlotTimeSeriesC() const
+{
+ return _isPlotTimeSeriesC;
+}
+
+bool ArgumentParser::GetIsPlotTimeSeriesD() const
+{
+ return _isPlotTimeSeriesD;
+}
+
+bool ArgumentParser::GetIsOneDimensional() const
+{
+ return _isOneDimensional;
+}
+
bool ArgumentParser::GetIsIndividualFD() const
{
return _isIndividualFD;
}
+polygon_2d ArgumentParser::GetAreaIndividualFD() const
+{
+ return _areaIndividualFD;
+}
+
bool ArgumentParser::GetIsGetProfile() const
{
return _isGetProfile;
@@ -682,7 +792,12 @@ vector<int> ArgumentParser::GetAreaIDforMethodD() const
MeasurementArea* ArgumentParser::GetMeasurementArea(int id)
{
if (_measurementAreas.count(id) == 0)
- return NULL;
+ {
+ Log->Write("ERROR:\t measurement id [%d] not found.",id);
+ Log->Write(" \t check your configuration files");
+ exit(EXIT_FAILURE);
+ //return NULL;
+ }
return _measurementAreas[id];
}
diff --git a/general/ArgumentParser.h b/general/ArgumentParser.h
index 03e230dc52c05d534133e68184ec9b9d0fc21070..19cc875fb10263efb3f0c0e18034266ff36dcb74 100644
--- a/general/ArgumentParser.h
+++ b/general/ArgumentParser.h
@@ -57,7 +57,7 @@ private:
FileFormat _fileFormat;
std::vector<std::string> _trajectoriesFiles;
- char _vComponent;
+ std::string _vComponent;
bool _isMethodA;
bool _isMethodB;
bool _isMethodC;
@@ -66,7 +66,13 @@ private:
double _cutRadius;
int _circleEdges;
bool _isOutputGraph;
+ bool _isPlotGraph;
+ bool _isPlotTimeSeriesA;
+ bool _isPlotTimeSeriesC;
+ bool _isPlotTimeSeriesD;
+ bool _isOneDimensional;
bool _isIndividualFD;
+ polygon_2d _areaIndividualFD;
bool _isGetProfile;
double _steadyStart;
double _steadyEnd;
@@ -103,7 +109,7 @@ public:
double GetLineEndX() const;
double GetLineEndY() const;
- char GetVComponent() const;
+ std::string GetVComponent() const;
int GetDelatT_Vins() const;
int GetTimeIntervalA() const;
bool GetIsMethodA() const;
@@ -118,7 +124,13 @@ public:
double GetCutRadius() const;
int GetCircleEdges() const;
bool GetIsOutputGraph() const;
+ bool GetIsPlotGraph() const;
+ bool GetIsPlotTimeSeriesA() const;
+ bool GetIsPlotTimeSeriesC() const;
+ bool GetIsPlotTimeSeriesD() const;
+ bool GetIsOneDimensional() const;
bool GetIsIndividualFD() const;
+ polygon_2d GetAreaIndividualFD() const;
double GetSteadyStart() const;
double GetSteadyEnd() const;
bool GetIsGetProfile() const;
diff --git a/geometry/Building.cpp b/geometry/Building.cpp
index 89400ce07315130e0c4cc5f877b8d4768b913703..2a9b44080db7d772d8b3f663bf8c7c8b5512ad27 100644
--- a/geometry/Building.cpp
+++ b/geometry/Building.cpp
@@ -892,7 +892,7 @@ bool Building::IsVisible(const Point& p1, const Point& p2, const std::vector<Sub
{
for(auto&& sub: subrooms)
{
- if(sub and sub->IsVisible(p1,p2,considerHlines)==false) return false;
+ if(sub && sub->IsVisible(p1,p2,considerHlines)==false) return false;
}
}
diff --git a/geometry/Line.cpp b/geometry/Line.cpp
index 1df68720eaedde5e6bcd6df33c7a72aa3ae16cb2..823e6d9af9cce9ac9bc14da07de2d3944361e453 100644
--- a/geometry/Line.cpp
+++ b/geometry/Line.cpp
@@ -268,13 +268,13 @@ bool Line::Overlapp(const Line& l) const
if(fabs(vecAB.Determinant(vecDC))<J_EPS)
{
- if( IsInLineSegment(l.GetPoint1()) and not HasEndPoint(l.GetPoint1()))
+ if( IsInLineSegment(l.GetPoint1()) && !HasEndPoint(l.GetPoint1()))
{
//Log->Write("ERROR: 1. Overlapping walls %s and %s ", toString().c_str(),l.toString().c_str());
return true;
}
- if( IsInLineSegment(l.GetPoint2()) and not HasEndPoint(l.GetPoint2()))
+ if( IsInLineSegment(l.GetPoint2()) && !HasEndPoint(l.GetPoint2()))
{
//Log->Write("ERROR: 2. Overlapping walls %s and %s ", toString().c_str(),l.toString().c_str());
return true;
diff --git a/geometry/SubRoom.cpp b/geometry/SubRoom.cpp
index 0b289cf90d84cca7ce0fb52854c435adef4edc53..37c9e02cc2aba3a7b295196bd2a81f76def56a8e 100644
--- a/geometry/SubRoom.cpp
+++ b/geometry/SubRoom.cpp
@@ -191,7 +191,7 @@ void SubRoom::AddTransition(Transition* line)
void SubRoom::AddNeighbor(SubRoom* sub)
{
- if(sub and (IsElementInVector(_neighbors, sub)==false))
+ if(sub && (IsElementInVector(_neighbors, sub)==false))
{
_neighbors.push_back(sub);
}
@@ -519,7 +519,7 @@ bool SubRoom::SanityCheck()
exit(EXIT_FAILURE);
//return false;
}
- connected=connected or w1.ShareCommonPointWith(w2);
+ connected=connected || w1.ShareCommonPointWith(w2);
}
//overlapping with lines
for(auto&& hline: _hlines)
@@ -540,7 +540,7 @@ bool SubRoom::SanityCheck()
exit(EXIT_FAILURE);
//return false;
}
- connected=connected or w1.ShareCommonPointWith(*c);
+ connected=connected || w1.ShareCommonPointWith(*c);
}
//overlaping with transitions
for(auto&& t: _transitions)
@@ -551,10 +551,10 @@ bool SubRoom::SanityCheck()
exit(EXIT_FAILURE);
//return false;
}
- connected=connected or w1.ShareCommonPointWith(*t);
+ connected=connected || w1.ShareCommonPointWith(*t);
}
- if(not connected)
+ if(!connected)
{
Log->Write("ERROR: loose wall found %s in Room/Subroom %d/%d",w1.toString().c_str(),_roomID,_id);
exit(EXIT_FAILURE);
diff --git a/methods/Method_A.cpp b/methods/Method_A.cpp
index bbfc37309b7e829a26e2dd07bc048c6e07f0eb81..630ffb5632706038aa87a11b5f0c91ac96b9809b 100644
--- a/methods/Method_A.cpp
+++ b/methods/Method_A.cpp
@@ -71,9 +71,11 @@ bool Method_A::Process (const PedData& peddata,const string& scriptsLocation)
_passLine[i] = false;
}
Log->Write("------------------------Analyzing with Method A-----------------------------");
- for(int frameNr = 0; frameNr < peddata.GetNumFrames(); frameNr++ )
+ //for(int frameNr = 0; frameNr < peddata.GetNumFrames(); frameNr++ )
+ for(std::map<int , std::vector<int> >::iterator ite=_peds_t.begin();ite!=_peds_t.end();ite++)
{
- int frid = frameNr + peddata.GetMinFrame();
+ int frameNr = ite->first;
+ int frid = frameNr + peddata.GetMinFrame();
if(!(frid%100))
{
Log->Write("frame ID = %d",frid);
@@ -100,11 +102,15 @@ void Method_A::WriteFile_N_t(string data)
if(file.is_open())
{
file<<data;
+ file.close();
string METHOD_A_LOCATION =_projectRootDir+"./Output/Fundamental_Diagram/FlowVelocity/";
string file_N_t ="Flow_NT_"+_trajName+"_id_"+_measureAreaId+".dat";
- string parameters_N_t="python "+_scriptsLocation+"/_Plot_N_t.py -p \""+ METHOD_A_LOCATION + "\" -n "+file_N_t;
- system(parameters_N_t.c_str());
- Log->Write("INFO:\tPlotting N-t diagram!");
+ if(_plotTimeSeries)
+ {
+ string parameters_N_t="python \""+_scriptsLocation+"/_Plot_N_t.py\" -p \""+ METHOD_A_LOCATION + "\" -n "+file_N_t;
+ int res = system(parameters_N_t.c_str());
+ Log->Write("INFO:\tPlotting N-t diagram! Status: %d", res);
+ }
}
else
{
@@ -213,3 +219,8 @@ void Method_A::SetTimeInterval(const int& deltaT)
{
_deltaT = deltaT;
}
+
+void Method_A::SetPlotTimeSeries(bool plotTimeseries)
+{
+ _plotTimeSeries = plotTimeseries;
+}
diff --git a/methods/Method_A.h b/methods/Method_A.h
index 62ea68bf05137eaf895f2c1693161104992c1f4f..4913aeef392e4151030f9c4367310d47bc803fb9 100644
--- a/methods/Method_A.h
+++ b/methods/Method_A.h
@@ -48,6 +48,7 @@ public:
void SetMeasurementArea (MeasurementArea_L* area);
void SetTimeInterval(const int& deltaT);
bool Process (const PedData& peddata,const std::string& scriptsLocation);
+ void SetPlotTimeSeries(bool plotTimeseries);
private:
std::string _trajName;
@@ -69,6 +70,8 @@ private:
int _classicFlow; // the number of pedestrians pass a line in a certain time
double _vDeltaT; // define this is to measure cumulative velocity each pedestrian pass a measure line each time step to calculate the <v>delat T=sum<vi>/N
int _deltaT;
+
+ bool _plotTimeSeries;
/**
* Calculate the Flow rate during a certain time interval DeltaT and the mean velocity passing a line.
* Note: here the time interval in calculating the flow rate is modified.
@@ -90,4 +93,4 @@ private:
void GetAccumFlowVelocity(int frame, const std::vector<int>& ids, const std::vector<double>& VInFrame);
};
-#endif /* METHOD_A_H_ */
\ No newline at end of file
+#endif /* METHOD_A_H_ */
diff --git a/methods/Method_B.cpp b/methods/Method_B.cpp
index 8b5c126c14a3c2cd19551c51e929d5926ad79b36..1e3fa6b0da3d4af172ce6617e1be6f9c94252e6d 100644
--- a/methods/Method_B.cpp
+++ b/methods/Method_B.cpp
@@ -132,7 +132,7 @@ void Method_B::GetFundamentalTinTout(double *DensityPerFrame,double LengthMeasur
fprintf(fFD_TinTout,"#person Index\t density_i(m^(-2))\t velocity_i(m/s)\n");
for(int i=0; i<_NumPeds; i++)
{
- double velocity_temp=_fps*CMtoM*LengthMeasurementarea/(_tOut[i]-_tIn[i]);
+ double velocity_temp=_fps*LengthMeasurementarea/(_tOut[i]-_tIn[i]);
double density_temp=0;
for(int j=_tIn[i]; j<_tOut[i]; j++)
{
diff --git a/methods/Method_C.cpp b/methods/Method_C.cpp
index 4a4b13ffbbc81b360935500c356525d4b64a4a26..0521d37895741306fc2eddc72fab623605fd5aa6 100644
--- a/methods/Method_C.cpp
+++ b/methods/Method_C.cpp
@@ -55,9 +55,11 @@ bool Method_C::Process (const PedData& peddata)
_fps = peddata.GetFps();
OpenFileMethodC();
Log->Write("------------------------Analyzing with Method C-----------------------------");
- for(int frameNr = 0; frameNr < peddata.GetNumFrames(); frameNr++ )
+ for(std::map<int , std::vector<int> >::iterator ite=_peds_t.begin();ite!=_peds_t.end();ite++)
{
- int frid = frameNr + _minFrame;
+ int frameNr = ite->first;
+ int frid = frameNr + _minFrame;
+
if(!(frid%100))
{
Log->Write("frame ID = %d",frid);
diff --git a/methods/Method_C.h b/methods/Method_C.h
index b1a1eddf1989d8ce49a2ea857e7a0fc9fb016e0d..d2018d76242156ac35aa22107b56d4479d66d043 100644
--- a/methods/Method_C.h
+++ b/methods/Method_C.h
@@ -55,4 +55,4 @@ private:
};
-#endif /* METHOD_C_H_ */
\ No newline at end of file
+#endif /* METHOD_C_H_ */
diff --git a/methods/Method_D.cpp b/methods/Method_D.cpp
index 43add369cf2b17761b0906c6c3439a843fa94969..7a1b9439b8fd2dcfea91440aed82ee3679e55195 100644
--- a/methods/Method_D.cpp
+++ b/methods/Method_D.cpp
@@ -40,6 +40,7 @@ Method_D::Method_D()
{
_grid_size_X = 0.10;
_grid_size_Y = 0.10;
+ _fps=16;
_outputVoronoiCellData = false;
_getProfile = false;
_geoMinX = 0;
@@ -53,6 +54,8 @@ Method_D::Method_D()
_calcIndividualFD = false;
_fVoronoiRhoV = NULL;
_areaForMethod_D = NULL;
+ _plotVoronoiCellData=false;
+ _isOneDimensional=false;
}
Method_D::~Method_D()
@@ -62,10 +65,6 @@ Method_D::~Method_D()
bool Method_D::Process (const PedData& peddata,const std::string& scriptsLocation)
{
- /*if(false==IsPedInGeometry(peddata.GetNumFrames(), peddata.GetNumPeds(), peddata.GetXCor(), peddata.GetYCor(), peddata.GetFirstFrame(), peddata.GetLastFrame()))
- {
- return false;
- }*/
_scriptsLocation = scriptsLocation;
_peds_t = peddata.GetPedsFrame();
_trajName = peddata.GetTrajName();
@@ -79,10 +78,11 @@ bool Method_D::Process (const PedData& peddata,const std::string& scriptsLocatio
OpenFileIndividualFD();
}
Log->Write("------------------------Analyzing with Method D-----------------------------");
- for(int frameNr = 0; frameNr < peddata.GetNumFrames(); frameNr++ )
+ //for(int frameNr = 0; frameNr < peddata.GetNumFrames(); frameNr++ )
+ for(std::map<int , std::vector<int> >::iterator ite=_peds_t.begin();ite!=_peds_t.end();ite++)
{
- int frid = frameNr + minFrame;
-
+ int frameNr = ite->first;
+ int frid = frameNr + minFrame;
//padd the frameid with 0
std::ostringstream ss;
ss << std::setw(5) << std::setfill('0') << frid;
@@ -116,23 +116,49 @@ bool Method_D::Process (const PedData& peddata,const std::string& scriptsLocatio
//---------------------------------------------------------------------------------------------------------------
if(NumPeds>2)
{
- vector<polygon_2d> polygons = GetPolygons(XInFrame, YInFrame, VInFrame, IdInFrame);
- OutputVoronoiResults(polygons, str_frid, VInFrame);
- if(_calcIndividualFD)
- {
- // if(i>beginstationary&&i<endstationary)
- {
- GetIndividualFD(polygons,VInFrame, IdInFrame, _areaForMethod_D->_poly, str_frid);
- }
- }
- if(_getProfile)
- { // field analysis
- GetProfiles(str_frid, polygons, VInFrame);
- }
- if(_outputVoronoiCellData)
- { // output the Voronoi polygons of a frame
- OutputVoroGraph(str_frid, polygons, NumPeds, XInFrame, YInFrame,VInFrame);
- }
+ if(_isOneDimensional)
+ {
+ CalcVoronoiResults1D(XInFrame, VInFrame, IdInFrame, _areaForMethod_D->_poly,str_frid);
+ }
+ else
+ {
+ if(IsPointsOnOneLine(XInFrame, YInFrame))
+ {
+ if(fabs(XInFrame[1]-XInFrame[0])<dmin)
+ {
+ XInFrame[1]+= offset;
+ }
+ else
+ {
+ YInFrame[1]+= offset;
+ }
+ }
+ vector<polygon_2d> polygons = GetPolygons(XInFrame, YInFrame, VInFrame, IdInFrame);
+ if(!polygons.empty())
+ {
+ OutputVoronoiResults(polygons, str_frid, VInFrame);
+ if(_calcIndividualFD)
+ {
+ if(false==_isOneDimensional)
+ {
+ GetIndividualFD(polygons,VInFrame, IdInFrame, _areaIndividualFD, str_frid);
+ }
+ }
+ if(_getProfile)
+ { // field analysis
+ GetProfiles(str_frid, polygons, VInFrame);
+ }
+ if(_outputVoronoiCellData)
+ { // output the Voronoi polygons of a frame
+ OutputVoroGraph(str_frid, polygons, NumPeds, XInFrame, YInFrame,VInFrame);
+ }
+ }
+ else
+ {
+ Log->Write("Warning:\tVoronoi Diagrams are not obtained!\n");
+ return false;
+ }
+ }
}
else
{
@@ -157,7 +183,14 @@ bool Method_D::OpenFileMethodD()
}
else
{
- fprintf(_fVoronoiRhoV,"#framerate:\t%.2f\n\n#Frame \t Voronoi density(m^(-2))\t Voronoi velocity(m/s)\n",_fps);
+ if(_isOneDimensional)
+ {
+ fprintf(_fVoronoiRhoV,"#framerate:\t%.2f\n\n#Frame \t Voronoi density(m^(-1))\t Voronoi velocity(m/s)\n",_fps);
+ }
+ else
+ {
+ fprintf(_fVoronoiRhoV,"#framerate:\t%.2f\n\n#Frame \t Voronoi density(m^(-2))\t Voronoi velocity(m/s)\n",_fps);
+ }
return true;
}
}
@@ -172,7 +205,14 @@ bool Method_D::OpenFileIndividualFD()
}
else
{
- fprintf(_fIndividualFD,"#framerate:\t%.2f\n\n#Frame \t PedId \t Individual density(m^(-2))\t Individual velocity(m/s)\n",_fps);
+ if(_isOneDimensional)
+ {
+ fprintf(_fIndividualFD,"#framerate:\t%.2f\n\n#Frame \t PedId \t headway(m)\t Individual velocity(m/s)\n",_fps);
+ }
+ else
+ {
+ fprintf(_fIndividualFD,"#framerate:\t%.2f\n\n#Frame \t PedId \t Individual density(m^(-2))\t Individual velocity(m/s)\n",_fps);
+ }
return true;
}
}
@@ -381,8 +421,7 @@ void Method_D::OutputVoroGraph(const string & frameId, vector<polygon_2d>& polyg
exit(EXIT_FAILURE);
}
-
- string point=voronoiLocation+"/points"+_trajName+"_id_"+_measureAreaId+"_"+frameId+".dat";
+ /*string point=voronoiLocation+"/points"+_trajName+"_id_"+_measureAreaId+"_"+frameId+".dat";
ofstream points (point.c_str());
if( points.is_open())
{
@@ -395,18 +434,20 @@ void Method_D::OutputVoroGraph(const string & frameId, vector<polygon_2d>& polyg
{
Log->Write("ERROR:\tcannot create the file <%s>",point.c_str());
exit(EXIT_FAILURE);
+ }*/
+
+ if(_plotVoronoiCellData)
+ {
+ string parameters_rho="python "+_scriptsLocation+"/_Plot_cell_rho.py -f \""+ voronoiLocation + "\" -n "+ _trajName+"_id_"+_measureAreaId+"_"+frameId+
+ " -g "+_geometryFileName;
+ string parameters_v="python "+_scriptsLocation+"/_Plot_cell_v.py -f \""+ voronoiLocation + "\" -n "+ _trajName+"_id_"+_measureAreaId+"_"+frameId+
+ " -g "+_geometryFileName;
+ //Log->Write("INFO:\t%s",parameters_rho.c_str());
+ Log->Write("INFO:\tPlotting Voronoi Cell at the frame <%s>",frameId.c_str());
+ system(parameters_rho.c_str());
+ system(parameters_v.c_str());
}
- string parameters_rho="python "+_scriptsLocation+"/_Plot_cell_rho.py -f \""+ voronoiLocation + "\" -n "+ _trajName+"_id_"+_measureAreaId+"_"+frameId+
- " -x1 "+boost::lexical_cast<std::string>(_geoMinX*CMtoM)+" -x2 "+boost::lexical_cast<std::string>(_geoMaxX*CMtoM)+" -y1 "+
- boost::lexical_cast<std::string>(_geoMinY*CMtoM)+" -y2 "+boost::lexical_cast<std::string>(_geoMaxY*CMtoM);
- string parameters_v="python "+_scriptsLocation+"/_Plot_cell_v.py -f \""+ voronoiLocation + "\" -n "+ _trajName+"_id_"+_measureAreaId+"_"+frameId+
- " -x1 "+boost::lexical_cast<std::string>(_geoMinX*CMtoM)+" -x2 "+boost::lexical_cast<std::string>(_geoMaxX*CMtoM)+" -y1 "+
- boost::lexical_cast<std::string>(_geoMinY*CMtoM)+" -y2 "+boost::lexical_cast<std::string>(_geoMaxY*CMtoM);
- //Log->Write("INFO:\t%s",parameters_rho.c_str());
- Log->Write("INFO:\tPlotting Voronoi Cell at the frame <%s>",frameId.c_str());
- system(parameters_rho.c_str());
- system(parameters_v.c_str());
- points.close();
+ //points.close();
polys.close();
velo.close();
}
@@ -437,6 +478,11 @@ void Method_D::SetCalculateIndividualFD(bool individualFD)
_calcIndividualFD = individualFD;
}
+void Method_D::SetAreaIndividualFD(polygon_2d areaindividualFD)
+{
+ _areaIndividualFD = areaindividualFD;
+}
+
void Method_D::Setcutbycircle(double radius,int edges)
{
_cutByCircle=true;
@@ -457,6 +503,11 @@ void Method_D::SetGeometryBoundaries(double minX, double minY, double maxX, doub
_geoMaxY = maxY;
}
+void Method_D::SetGeometryFileName(const std::string& geometryFile)
+{
+ _geometryFileName=geometryFile;
+}
+
void Method_D::SetGridSize(double x, double y)
{
_grid_size_X = x;
@@ -473,11 +524,21 @@ void Method_D::SetOutputVoronoiCellData(bool outputCellData)
_outputVoronoiCellData = outputCellData;
}
+void Method_D::SetPlotVoronoiGraph(bool plotVoronoiGraph)
+{
+ _plotVoronoiCellData = plotVoronoiGraph;
+}
+
void Method_D::SetMeasurementArea (MeasurementArea_B* area)
{
_areaForMethod_D = area;
}
+void Method_D::SetDimensional (bool dimension)
+{
+ _isOneDimensional = dimension;
+}
+
void Method_D::ReducePrecision(polygon_2d& polygon)
{
for(auto&& point:polygon.outer())
@@ -500,3 +561,121 @@ bool Method_D::IsPedInGeometry(int frames, int peds, double **Xcor, double **Yco
}
return true;
}
+
+void Method_D::CalcVoronoiResults1D(vector<double>& XInFrame, vector<double>& VInFrame, vector<int>& IdInFrame, const polygon_2d & measureArea,const string& frid)
+{
+ vector<double> measurearea_x;
+ for(unsigned int i=0;i<measureArea.outer().size();i++)
+ {
+ measurearea_x.push_back(measureArea.outer()[i].x());
+ }
+ double left_boundary=*min_element(measurearea_x.begin(),measurearea_x.end());
+ double right_boundary=*max_element(measurearea_x.begin(),measurearea_x.end());
+
+ vector<double> voronoi_distance;
+ vector<double> Xtemp=XInFrame;
+ vector<double> dist;
+ vector<double> XLeftNeighbor;
+ vector<double> XLeftTemp;
+ vector<double> XRightNeighbor;
+ vector<double> XRightTemp;
+ sort(Xtemp.begin(),Xtemp.end());
+ dist.push_back(Xtemp[1]-Xtemp[0]);
+ XLeftTemp.push_back(2*Xtemp[0]-Xtemp[1]);
+ XRightTemp.push_back(Xtemp[1]);
+ for(unsigned int i=1;i<Xtemp.size()-1;i++)
+ {
+ dist.push_back((Xtemp[i+1]-Xtemp[i-1])/2.0);
+ XLeftTemp.push_back(Xtemp[i-1]);
+ XRightTemp.push_back(Xtemp[i+1]);
+ }
+ dist.push_back(Xtemp[Xtemp.size()-1]-Xtemp[Xtemp.size()-2]);
+ XLeftTemp.push_back(Xtemp[Xtemp.size()-2]);
+ XRightTemp.push_back(2*Xtemp[Xtemp.size()-1]-Xtemp[Xtemp.size()-2]);
+ for(unsigned int i=0;i<XInFrame.size();i++)
+ {
+ for(unsigned int j=0;j<Xtemp.size();j++)
+ {
+ if(fabs(XInFrame[i]-Xtemp[j])<1.0e-5)
+ {
+ voronoi_distance.push_back(dist[j]);
+ XLeftNeighbor.push_back(XLeftTemp[j]);
+ XRightNeighbor.push_back(XRightTemp[j]);
+ break;
+ }
+ }
+ }
+
+ double VoronoiDensity=0;
+ double VoronoiVelocity=0;
+ for(unsigned int i=0; i<XInFrame.size(); i++)
+ {
+ double ratio=getOverlapRatio((XInFrame[i]+XLeftNeighbor[i])/2.0, (XRightNeighbor[i]+XInFrame[i])/2.0,left_boundary,right_boundary);
+ VoronoiDensity+=ratio;
+ VoronoiVelocity+=(VInFrame[i]*voronoi_distance[i]*ratio*CMtoM);
+ if(_calcIndividualFD)
+ {
+ double headway=(XRightNeighbor[i] - XInFrame[i])*CMtoM;
+ fprintf(_fIndividualFD,"%s\t%d\t%.3f\t%.3f\n",frid.c_str(), IdInFrame[i], headway,VInFrame[i]);
+ }
+ }
+ VoronoiDensity/=((right_boundary-left_boundary)*CMtoM);
+ VoronoiVelocity/=((right_boundary-left_boundary)*CMtoM);
+ fprintf(_fVoronoiRhoV,"%s\t%.3f\t%.3f\n",frid.c_str(),VoronoiDensity, VoronoiVelocity);
+
+}
+
+double Method_D::getOverlapRatio(const double& left, const double& right, const double& measurearea_left, const double& measurearea_right)
+{
+ double OverlapRatio=0;
+ double PersonalSpace=right-left;
+ if(left > measurearea_left && right < measurearea_right) //case1
+ {
+ OverlapRatio=1;
+ }
+ else if(right > measurearea_left && right < measurearea_right && left < measurearea_left)
+ {
+ OverlapRatio=(right-measurearea_left)/PersonalSpace;
+ }
+ else if(left < measurearea_left && right > measurearea_right)
+ {
+ OverlapRatio=(measurearea_right - measurearea_left)/PersonalSpace;
+ }
+ else if(left > measurearea_left && left < measurearea_right && right > measurearea_right)
+ {
+ OverlapRatio=(measurearea_right-left)/PersonalSpace;
+ }
+ return OverlapRatio;
+}
+
+bool Method_D::IsPointsOnOneLine(vector<double>& XInFrame, vector<double>& YInFrame)
+{
+ double deltaX=XInFrame[1] - XInFrame[0];
+ bool isOnOneLine=true;
+ if(fabs(deltaX)<dmin)
+ {
+ for(unsigned int i=2; i<XInFrame.size();i++)
+ {
+ if(fabs(XInFrame[i] - XInFrame[0])> dmin)
+ {
+ isOnOneLine=false;
+ break;
+ }
+ }
+ }
+ else
+ {
+ double slope=(YInFrame[1] - YInFrame[0])/deltaX;
+ double intercept=YInFrame[0] - slope*XInFrame[0];
+ for(unsigned int i=2; i<XInFrame.size();i++)
+ {
+ double dist=fabs(slope*XInFrame[i] - YInFrame[i] + intercept)/sqrt(slope*slope +1);
+ if(dist > dmin)
+ {
+ isOnOneLine=false;
+ break;
+ }
+ }
+ }
+ return isOnOneLine;
+}
diff --git a/methods/Method_D.h b/methods/Method_D.h
index 95b05f95d8ab9add585e6640bd72be0d9e63c9a1..8bf9ea5971af71714eb88aa976dc1ebd482d027a 100644
--- a/methods/Method_D.h
+++ b/methods/Method_D.h
@@ -31,7 +31,7 @@
#include "PedData.h"
#include "../Analysis.h"
#include "VoronoiDiagram.h"
-
+#include <algorithm>
#ifdef __linux__
#include <sys/stat.h>
@@ -43,6 +43,9 @@
#include <direct.h>
#endif
+//handle more than two person are in one line
+#define dmin 200
+#define offset 200
class Method_D
@@ -52,13 +55,17 @@ public:
virtual ~Method_D();
bool Process (const PedData& peddata,const std::string& scriptsLocation);
void SetCalculateIndividualFD(bool individualFD);
+ void SetAreaIndividualFD(polygon_2d areaindividualFD);
void Setcutbycircle(double radius,int edges);
void SetGeometryPolygon(polygon_2d geometryPolygon);
+ void SetGeometryFileName(const std::string& geometryFile);
void SetGeometryBoundaries(double minX, double minY, double maxX, double maxY);
void SetGridSize(double x, double y);
void SetCalculateProfiles(bool calcProfile);
void SetOutputVoronoiCellData(bool outputCellData);
+ void SetPlotVoronoiGraph(bool plotVoronoiGraph);
void SetMeasurementArea (MeasurementArea_B* area);
+ void SetDimensional (bool dimension);
private:
std::map<int , std::vector<int> > _peds_t;
@@ -68,8 +75,11 @@ private:
std::string _projectRootDir;
std::string _scriptsLocation;
bool _calcIndividualFD;
+ polygon_2d _areaIndividualFD;
bool _getProfile;
bool _outputVoronoiCellData;
+ bool _plotVoronoiCellData;
+ bool _isOneDimensional;
bool _cutByCircle; //Adjust whether cut each original voronoi cell by a circle
double _cutRadius;
int _circleEdges;
@@ -85,6 +95,7 @@ private:
float _fps;
bool OpenFileMethodD();
bool OpenFileIndividualFD();
+ std::string _geometryFileName;
std::vector<polygon_2d> GetPolygons(std::vector<double>& XInFrame, std::vector<double>& YInFrame,
std::vector<double>& VInFrame, std::vector<int>& IdInFrame);
@@ -96,13 +107,15 @@ private:
void OutputVoroGraph(const std::string & frameId, std::vector<polygon_2d>& polygons, int numPedsInFrame,std::vector<double>& XInFrame,
std::vector<double>& YInFrame,const std::vector<double>& VInFrame);
void GetIndividualFD(const std::vector<polygon_2d>& polygon, const std::vector<double>& Velocity, const std::vector<int>& Id, const polygon_2d& measureArea, const std::string& frid);
-
/**
* Reduce the precision of the points to two digits
* @param polygon
*/
+ void CalcVoronoiResults1D(std::vector<double>& XInFrame, std::vector<double>& VInFrame, std::vector<int>& IdInFrame, const polygon_2d & measureArea, const std::string& frid);
void ReducePrecision(polygon_2d& polygon);
bool IsPedInGeometry(int frames, int peds, double **Xcor, double **Ycor, int *firstFrame, int *lastFrame); //check whether all the pedestrians are in the geometry
+ double getOverlapRatio(const double& left, const double& right, const double& measurearea_left, const double& measurearea_right);
+ bool IsPointsOnOneLine(std::vector<double>& XInFrame, std::vector<double>& YInFrame);
};
-#endif /* METHOD_D_H_ */
\ No newline at end of file
+#endif /* METHOD_D_H_ */
diff --git a/methods/PedData.cpp b/methods/PedData.cpp
index b6a9c0d80d1fff8bdb342e2a24c6b8fe1ef3a620..185ba88992c423ddd7a3028505349af90a45b03e 100644
--- a/methods/PedData.cpp
+++ b/methods/PedData.cpp
@@ -1,7 +1,7 @@
/**
* \file PedData.cpp
- * \date Oct 10, 2014
- * \version v0.7
+ * \date Feb 10, 2016
+ * \version v0.8
* \copyright <2009-2015> Forschungszentrum J��lich GmbH. All rights reserved.
*
* \section License
@@ -44,7 +44,7 @@ PedData::~PedData()
}
-bool PedData::ReadData(const string& projectRootDir, const string& path, const string& filename, const FileFormat& trajformat, int deltaF, char vComponent)
+bool PedData::ReadData(const string& projectRootDir, const string& path, const string& filename, const FileFormat& trajformat, int deltaF, std::string vComponent)
{
_minID = INT_MAX;
_minFrame = INT_MAX;
@@ -80,6 +80,7 @@ bool PedData::InitializeVariables(const string& filename)
{
vector<double> xs;
vector<double> ys;
+ vector<string> vcmp; // the direction identification for velocity calculation
vector<int> _IdsTXT; // the Id data from txt format trajectory data
vector<int> _FramesTXT; // the Frame data from txt format trajectory data
//string fullTrajectoriesPathName= _projectRootDir+"./"+_trajName;
@@ -115,7 +116,7 @@ bool PedData::InitializeVariables(const string& filename)
std::vector<std::string> strs;
boost::split(strs, line , boost::is_any_of("\t "),boost::token_compress_on);
- if(strs.size() <4)
+ if(strs.size() < 4)
{
Log->Write("ERROR:\t There is an error in the file at line %d", lineNr);
return false;
@@ -125,6 +126,18 @@ bool PedData::InitializeVariables(const string& filename)
_FramesTXT.push_back(atoi(strs[1].c_str()));
xs.push_back(atof(strs[2].c_str()));
ys.push_back(atof(strs[3].c_str()));
+ if(_vComponent=="F")
+ {
+ if(strs.size() >= 6)
+ {
+ vcmp.push_back(strs[5].c_str());
+ }
+ else
+ {
+ Log->Write("ERROR:\t There is no indicator for velocity component in trajectory file or ini file!!");
+ return false;
+ }
+ }
}
lineNr++;
}
@@ -137,16 +150,10 @@ bool PedData::InitializeVariables(const string& filename)
//Total number of frames
_numFrames = *max_element(_FramesTXT.begin(),_FramesTXT.end()) - _minFrame+1;
+
//Total number of agents
_numPeds = *max_element(_IdsTXT.begin(),_IdsTXT.end()) - _minID+1;
- vector<int> Ids_temp=_IdsTXT;
- sort (Ids_temp.begin(), Ids_temp.end());
- Ids_temp.erase(unique(Ids_temp.begin(), Ids_temp.end()), Ids_temp.end());
- if((unsigned)_numPeds!=Ids_temp.size())
- {
- Log->Write("Error:\tThe index of ped ID is not continuous. Please modify the trajectory file!");
- return false;
- }
+
CreateGlobalVariables(_numPeds, _numFrames);
std::vector<int> firstFrameIndex; //The first frame index of each pedestrian
@@ -165,10 +172,18 @@ bool PedData::InitializeVariables(const string& filename)
lastFrameIndex.push_back(firstFrameIndex[i] - 1);
}
lastFrameIndex.push_back(_IdsTXT.size() - 1);
+
for (unsigned int i = 0; i < firstFrameIndex.size(); i++)
{
_firstFrame[i] = _FramesTXT[firstFrameIndex[i]] - _minFrame;
_lastFrame[i] = _FramesTXT[lastFrameIndex[i]] - _minFrame;
+ int actual_totalframe=lastFrameIndex[i]-firstFrameIndex[i]+1;
+ int expect_totalframe=_lastFrame[i]-_firstFrame[i]+1;
+ if(actual_totalframe != expect_totalframe)
+ {
+ Log->Write("Error:\tThe trajectory of ped with ID <%d> is not continuous. Please modify the trajectory file!",_IdsTXT[firstFrameIndex[i]]);
+ return false;
+ }
}
for(unsigned int i = 0; i < _IdsTXT.size(); i++)
@@ -179,6 +194,14 @@ bool PedData::InitializeVariables(const string& filename)
double y = ys[i]*M2CM;
_xCor[ID][frm] = x;
_yCor[ID][frm] = y;
+ if(_vComponent == "F")
+ {
+ _vComp[ID][frm] = vcmp[i];
+ }
+ else
+ {
+ _vComp[ID][frm] = _vComponent;
+ }
}
//save the data for each frame
@@ -187,6 +210,7 @@ bool PedData::InitializeVariables(const string& filename)
int id = _IdsTXT[i]-_minID;
int t =_FramesTXT[i]- _minFrame;
_peds_t[t].push_back(id);
+
}
return true;
@@ -204,15 +228,6 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
return false;
}
- //counting the number of frames
- int frames = 0;
- for(TiXmlElement* xFrame = xRootNode->FirstChildElement("frame"); xFrame;
- xFrame = xFrame->NextSiblingElement("frame")) {
- frames++;
- }
- _numFrames = frames;
- Log->Write("INFO:\tnum Frames = %d",_numFrames);
-
//Number of agents
TiXmlNode* xHeader = xRootNode->FirstChild("header"); // header
@@ -227,7 +242,6 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
Log->Write("INFO:\tFrame rate fps: <%.2f>", _fps);
}
- CreateGlobalVariables(_numPeds, _numFrames);
//processing the frames node
TiXmlNode* xFramesNode = xRootNode->FirstChild("frame");
@@ -237,6 +251,7 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
}
// obtaining the minimum id and minimum frame
+ int maxFrame=0;
for(TiXmlElement* xFrame = xRootNode->FirstChildElement("frame"); xFrame;
xFrame = xFrame->NextSiblingElement("frame"))
{
@@ -245,6 +260,10 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
{
_minFrame = frm;
}
+ if(frm>maxFrame)
+ {
+ maxFrame=frm;
+ }
for(TiXmlElement* xAgent = xFrame->FirstChildElement("agent"); xAgent;
xAgent = xAgent->NextSiblingElement("agent"))
{
@@ -255,10 +274,22 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
}
}
}
- int frameNr=0;
+
+ //counting the number of frames
+ _numFrames = maxFrame-_minFrame+1;
+ Log->Write("INFO:\tnum Frames = %d",_numFrames);
+
+ CreateGlobalVariables(_numPeds, _numFrames);
+
+ vector<int> totalframes;
+ for (int i = 0; i <_numPeds; i++)
+ {
+ totalframes.push_back(0);
+ }
+ //int frameNr=0;
for(TiXmlElement* xFrame = xRootNode->FirstChildElement("frame"); xFrame;
xFrame = xFrame->NextSiblingElement("frame")) {
-
+ int frameNr = atoi(xFrame->Attribute("ID")) - _minFrame;
//todo: can be parallelized with OpenMP
for(TiXmlElement* xAgent = xFrame->FirstChildElement("agent"); xAgent;
xAgent = xAgent->NextSiblingElement("agent")) {
@@ -268,9 +299,27 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
double y= atof(xAgent->Attribute("y"));
int ID= atoi(xAgent->Attribute("ID"))-_minID;
+
_peds_t[frameNr].push_back(ID);
_xCor[ID][frameNr] = x*M2CM;
_yCor[ID][frameNr] = y*M2CM;
+ if(_vComponent == "F")
+ {
+ if(xAgent->Attribute("vc"))
+ {
+ _vComp[ID][frameNr] = *string(xAgent->Attribute("vc")).c_str();
+ }
+ else
+ {
+ Log->Write("ERROR:\t There is no indicator for velocity component in trajectory file or ini file!!");
+ return false;
+ }
+ }
+ else
+ {
+ _vComp[ID][frameNr] = _vComponent;
+ }
+
if(frameNr < _firstFrame[ID])
{
_firstFrame[ID] = frameNr;
@@ -279,9 +328,22 @@ bool PedData::InitializeVariables(TiXmlElement* xRootNode)
{
_lastFrame[ID] = frameNr;
}
+ totalframes[ID] +=1;
}
- frameNr++;
+ //frameNr++;
+ }
+
+ for(int id = 0; id<_numPeds; id++)
+ {
+ int actual_totalframe= totalframes[id];
+ int expect_totalframe=_lastFrame[id]-_firstFrame[id]+1;
+ if(actual_totalframe != expect_totalframe)
+ {
+ Log->Write("Error:\tThe trajectory of ped with ID <%d> is not continuous. Please modify the trajectory file!",id+_minID);
+ return false;
+ }
}
+
return true;
}
@@ -333,10 +395,10 @@ vector<int> PedData::GetIdInFrame(const vector<int>& ids) const
double PedData::GetInstantaneousVelocity(int Tnow,int Tpast, int Tfuture, int ID, int *Tfirst, int *Tlast, double **Xcor, double **Ycor) const
{
-
+ std::string vcmp = _vComp[ID][Tnow];
double v=0.0;
-
- if(_vComponent == 'X')
+ //check the component used in the calculation of velocity
+ if(vcmp == "X" || vcmp == "X+"|| vcmp == "X-")
{
if((Tpast >=Tfirst[ID])&&(Tfuture <= Tlast[ID]))
{
@@ -350,8 +412,10 @@ double PedData::GetInstantaneousVelocity(int Tnow,int Tpast, int Tfuture, int ID
{
v = _fps*CMtoM*(Xcor[ID][Tnow] - Xcor[ID][Tpast])/( _deltaF); //one dimensional velocity
}
+ if((vcmp=="X+"&& v<0)||(vcmp=="X-"&& v>0)) //no moveback
+ v=0;
}
- if(_vComponent == 'Y')
+ else if(vcmp == "Y" || vcmp == "Y+"|| vcmp == "Y-")
{
if((Tpast >=Tfirst[ID])&&(Tfuture <= Tlast[ID]))
{
@@ -365,8 +429,10 @@ double PedData::GetInstantaneousVelocity(int Tnow,int Tpast, int Tfuture, int ID
{
v = _fps*CMtoM*(Ycor[ID][Tnow] - Ycor[ID][Tpast])/( _deltaF); //one dimensional velocity
}
+ if((vcmp=="Y+"&& v<0)||(vcmp=="Y-"&& v>0)) //no moveback
+ v=0;
}
- if(_vComponent == 'B')
+ else if(vcmp == "B")
{
if((Tpast >=Tfirst[ID])&&(Tfuture <= Tlast[ID]))
{
@@ -389,9 +455,11 @@ void PedData::CreateGlobalVariables(int numPeds, int numFrames)
{
_xCor = new double* [numPeds];
_yCor = new double* [numPeds];
+ _vComp = new string* [numPeds];
for (int i=0; i<numPeds; i++) {
_xCor[i] = new double [numFrames];
_yCor[i] = new double [numFrames];
+ _vComp[i] =new string [numFrames];
}
_firstFrame = new int[numPeds]; // Record the first frame of each pedestrian
_lastFrame = new int[numPeds]; // Record the last frame of each pedestrian
@@ -400,6 +468,7 @@ void PedData::CreateGlobalVariables(int numPeds, int numFrames)
for (int j = 0; j < numFrames; j++) {
_xCor[i][j] = 0;
_yCor[i][j] = 0;
+ _vComp[i][j] ="B";
}
_firstFrame[i] = INT_MAX;
_lastFrame[i] = INT_MIN;
diff --git a/methods/PedData.h b/methods/PedData.h
index 74d9341668a45eda4d11e54a668c1ca0156549ba..298181a119545eb3748ebb6a33b8f53113275f85 100644
--- a/methods/PedData.h
+++ b/methods/PedData.h
@@ -65,7 +65,7 @@ public:
std::vector<double> GetXInFrame(int frame, const std::vector<int>& ids) const;
std::vector<double> GetYInFrame(int frame, const std::vector<int>& ids) const;
std::vector<double> GetVInFrame(int frame, const std::vector<int>& ids) const;
- bool ReadData(const std::string& projectRootDir, const std::string& path, const std::string& filename, const FileFormat& _trajformat, int deltaF, char vComponent);
+ bool ReadData(const std::string& projectRootDir, const std::string& path, const std::string& filename, const FileFormat& _trajformat, int deltaF, std::string vComponent);
private:
@@ -87,13 +87,14 @@ private:
std::map<int , std::vector<int>> _peds_t;
int _deltaF=5;
- char _vComponent='X';
+ std::string _vComponent="B";
int *_firstFrame=NULL; // Record the first frame of each pedestrian
int *_lastFrame=NULL; // Record the last frame of each pedestrian
double **_xCor=NULL;
double **_yCor=NULL;
+ std::string **_vComp=NULL;
};
-#endif /* PEDDATA_H_ */
\ No newline at end of file
+#endif /* PEDDATA_H_ */
diff --git a/methods/VoronoiDiagram.cpp b/methods/VoronoiDiagram.cpp
index e8ed8314512741455445f44a5e800a97498b6781..89ca23fdd94ef1b0b752cfea46f84328a0e1d457 100644
--- a/methods/VoronoiDiagram.cpp
+++ b/methods/VoronoiDiagram.cpp
@@ -49,19 +49,19 @@ VoronoiDiagram::~VoronoiDiagram()
vector<polygon_2d> VoronoiDiagram::getVoronoiPolygons(vector<double>& XInFrame, vector<double>& YInFrame,
vector<double>& VInFrame, vector<int>& IdInFrame, const double Bound_Max)
{
- int numPedsInFrame = IdInFrame.size();
- int XInFrame_temp[numPedsInFrame];
- int YInFrame_temp[numPedsInFrame];
- double VInFrame_temp[numPedsInFrame];
- int IdInFrame_temp[numPedsInFrame];
+ const int numPedsInFrame = IdInFrame.size();
+ vector<int> XInFrame_temp;
+ vector<int> YInFrame_temp;
+ vector<double> VInFrame_temp;
+ vector<int> IdInFrame_temp;
for (int i = 0; i < numPedsInFrame; i++)
{
points.push_back(point_type2(round(XInFrame[i]), round(YInFrame[i])));
- XInFrame_temp[i] = round(XInFrame[i]);
- YInFrame_temp[i] = round(YInFrame[i]);
- VInFrame_temp[i] = VInFrame[i];
- IdInFrame_temp[i] = IdInFrame[i];
+ XInFrame_temp.push_back(round(XInFrame[i]));
+ YInFrame_temp.push_back(round(YInFrame[i]));
+ VInFrame_temp.push_back(VInFrame[i]);
+ IdInFrame_temp.push_back(IdInFrame[i]);
}
VD voronoidiagram;
@@ -83,7 +83,7 @@ vector<polygon_2d> VoronoiDiagram::getVoronoiPolygons(vector<double>& XInFrame,
for (voronoi_diagram<double>::const_cell_iterator it = voronoidiagram.cells().begin();
it != voronoidiagram.cells().end(); ++it)
{
- polygon_2d poly;
+ polygon_2d poly;
vector<point_type2> polypts;
point_type2 pt_s;
point_type2 pt_e;
@@ -277,7 +277,7 @@ bool VoronoiDiagram::point_inside_triangle(const point_type2& pt, const point_ty
vector<point_type2> VoronoiDiagram::add_bounding_points(const point_type2& pt1, const point_type2& pt2, const point_type2& pt,
double minX, double minY, double maxX, double maxY)
{
- vector<point_type2> bounding_vertex;
+ vector<point_type2> bounding_vertex;
if (fabs(pt1.x() - pt2.x()) > J_EPS && fabs(pt1.y() - pt2.y()) > J_EPS)
{
if ((fabs(pt1.y() - maxY) < J_EPS && fabs(pt2.x() - maxX) < J_EPS)
@@ -359,7 +359,7 @@ vector<point_type2> VoronoiDiagram::add_bounding_points(const point_type2& pt1,
else
{
double tempx = pt1.x()
- + (pt2.x() - pt1.x()) * (pt.y() - pt1.y()) / (pt2.y() - pt1.y());
+ + (pt2.x() - pt1.x()) * (pt.y() - pt1.y()) / (pt2.y() - pt1.y());
if (tempx < pt.x())
{
bounding_vertex.push_back(point_type2(minX, maxY));
@@ -391,7 +391,7 @@ vector<point_type2> VoronoiDiagram::add_bounding_points(const point_type2& pt1,
else
{
double tempy = pt1.y()
- + (pt2.y() - pt1.y()) * (pt.x() - pt1.x()) / (pt2.x() - pt1.x());
+ + (pt2.y() - pt1.y()) * (pt.x() - pt1.x()) / (pt2.x() - pt1.x());
if (tempy < pt.y())
{
bounding_vertex.push_back(point_type2(minX, minY));
@@ -436,7 +436,7 @@ point_type2 VoronoiDiagram::getIntersectionPoint(const point_2d& pt0, const poin
std::vector<polygon_2d> VoronoiDiagram::cutPolygonsWithGeometry(const std::vector<polygon_2d>& polygon,
const polygon_2d& Geometry, const vector<double>& xs, const vector<double>& ys)
{
- vector<polygon_2d> intersetionpolygons;
+ vector<polygon_2d> intersetionpolygons;
int temp = 0;
for (const auto & polygon_iterator:polygon)
{
@@ -464,7 +464,7 @@ std::vector<polygon_2d> VoronoiDiagram::cutPolygonsWithGeometry(const std::vecto
std::vector<polygon_2d> VoronoiDiagram::cutPolygonsWithCircle(const std::vector<polygon_2d>& polygon,
const vector<double>& xs, const vector<double>& ys, double radius, int edges)
{
- std::vector<polygon_2d> intersetionpolygons;
+ std::vector<polygon_2d> intersetionpolygons;
int temp = 0;
for (const auto & polygon_iterator:polygon)
{
@@ -482,7 +482,7 @@ std::vector<polygon_2d> VoronoiDiagram::cutPolygonsWithCircle(const std::vector<
intersection(circle, polygon_iterator, v);
BOOST_FOREACH(auto & it, v)
{
- if (within(point_2d(xs[temp], ys[temp]), it))
+ if (within(point_2d(xs[temp], ys[temp]), it))
{
//check and remove duplicates
//dispatch::unique (it);
@@ -496,4 +496,3 @@ std::vector<polygon_2d> VoronoiDiagram::cutPolygonsWithCircle(const std::vector<
}
return intersetionpolygons;
}
-
diff --git a/scripts/SteadyState.py b/scripts/SteadyState.py
index 681d3c8db83a2ae3e4f535833a8abd9fbaa51031..7a56f6a4947e15e40a625e4c594dc12ebd177c28 100644
--- a/scripts/SteadyState.py
+++ b/scripts/SteadyState.py
@@ -1,11 +1,17 @@
+import os
+from sys import exit
import argparse
import scipy.stats
from numpy import *
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches
from scipy.stats.stats import pearsonr
+import itertools
-frame = 16
+VERSION = 0.8
+script_dir = os.path.dirname(os.path.realpath(__file__))
+demos_dir = os.path.join(os.path.dirname(os.path.abspath(script_dir)), "demos")
+demo_file = os.path.join(demos_dir, "SteadyState", "rho_v_Voronoi_traj_AO_b240.txt_id_1.dat")
alpha = 0.99
gamma = 0.99
@@ -16,413 +22,405 @@ dx = 0.01
NN = s_max
KK = 500
im = 3.2
-d = 2*im/KK
-xi = arange(-im, im+0.0001, d)
-ia = ib = 0
-dnorm = 0
+d = 2 * im / KK
+xi = arange(-im, im + 0.0001, d)
+
def F(x):
+ """
+ :param x:
+ :return:
+ """
if x > q:
return 1
else:
return -1
-
-def func_b(i, k):
- return d * exp(-(xi[i]-xi[k]*acf) * (xi[i]-xi[k]*acf) / (2*(1-acf*acf))) / sqrt(2*math.pi*(1-acf*acf))
+
+def func_b(i, k, acf):
+ """
+ :param i:
+ :param k:
+ :param acf:
+ :return:
+ """
+ return d * exp(-(xi[i] - xi[k] * acf) * (xi[i] - xi[k] * acf) / (2 * (1 - acf * acf))) / sqrt(
+ 2 * math.pi * (1 - acf * acf))
+
+
+# - Flow
+# - plot frame x axis
+# - Unities /
def getParserArgs():
- parser = argparse.ArgumentParser(description='Combine French data to one file')
- parser.add_argument("-v", "--version", help='JuPedSim Version 0.7 Detection of Steady State')
- parser.add_argument("-f", "--input_file", default="E:/GitLab/jpsreport/demos/SteadyState/rho_v_Voronoi_traj_AO_b240.txt_id_1.dat", help='give the location and the name of the input file')
- parser.add_argument("-rs", "--reference_rho_start", type=int, default=240, help='give the start frame of the reference process in density')
- parser.add_argument("-re", "--reference_rho_end", type=int, default=640, help='give the end frame of the reference process in density')
- parser.add_argument("-vs", "--reference_v_start", type=int, default=240, help='give the start frame of the reference process in speed')
- parser.add_argument("-ve", "--reference_v_end", type=int, default=640, help='give the end frame of the reference process in speed')
- parser.add_argument("-p", "--plotfigs", default="yes", help='choose to plot the figures or not')
+ """
+
+ :return:
+ """
+ parser = argparse.ArgumentParser(
+ description='Detection of Steady State. JuPedSim v%0.1f.' % VERSION)
+ parser.add_argument("-a", "--automatic", action='store_true',
+ help="use calculated reference start and end frames. If given -rs and -re are ignored")
+ parser.add_argument("-f", "--input_file", default=demo_file,
+ help='Full name of the input file (default %s). File format is | frame | rho | velocity |' % demo_file)
+ parser.add_argument("-c", "--columns", nargs='+', type=int, default=[0, 1, 2],
+ help="columns to read from file (default 0, 1, 2)")
+ parser.add_argument("-xl", "--xlabel", type=str, default="t /s",
+ help="xlabel")
+ parser.add_argument("-yl", "--ylabel", nargs='+', type=str, default=["\rho /m/s", "v m/s"],
+ help="ylabel")
+ parser.add_argument("-rs", "--reference_start", nargs='+', type=int, default=[240, 240],
+ help='Start frame of the reference process in density (default 240)')
+ parser.add_argument("-re", "--reference_end", nargs='+', type=int, default=[640, 640],
+ help='End frame of the reference process in density (default 640)')
+ parser.add_argument("-p", "--plotfigs", default="yes",
+ help='Plot figures (default yes)')
+ parser.add_argument("-r", "--fps", type=int, default=16,
+ help='Frame per second (default 16)')
args = parser.parse_args()
return args
-if __name__ == '__main__':
- rho_max = 8.0
- args = getParserArgs()
- input_file = args.input_file
- ref_rho_start = args.reference_rho_start
- ref_rho_end = args.reference_rho_end
- ref_v_start = args.reference_v_start
- ref_v_end = args.reference_v_end
- plotfigs = args.plotfigs
-
-# input data
-data = loadtxt('%s'%(input_file))
-data = data[ data[:,1] != 0 ]
-data[:,0] = data[:,0] - data[0,0]
-
-# filepath and filename
-reverse = input_file[::-1]
-count = -1
-for i in reverse:
- count += 1
- if i == '/':
- break
-filepath = input_file[:(len(input_file)-count)]
-filename = input_file[-count:]
-print('file path = ', filepath)
-print('file name = ', filename)
-
-# calculate ref_mean, ref_std, ref_acf for rho
-ref_rho = data
-ref_rho = ref_rho[ ref_rho[:,0] >= ref_rho_start ]
-ref_rho = ref_rho[ ref_rho[:,0] <= ref_rho_end ]
-ref_rho = ref_rho[:,1]
-ref_rho_mean = mean(ref_rho)
-ref_rho_std = std(ref_rho)
-ref_rho_a = ref_rho[:-1]
-ref_rho_b = ref_rho[1:]
-ref_rho_acf = pearsonr(ref_rho_a, ref_rho_b)
-ref_rho_acf = ref_rho_acf[0]
-
-# calculate ref_mean, ref_std, ref_acf for v
-ref_v = data
-ref_v = ref_v[ ref_v[:,0] >= ref_v_start ]
-ref_v = ref_v[ ref_v[:,0] <= ref_v_end ]
-ref_v = ref_v[:,2]
-ref_v_mean = mean(ref_v)
-ref_v_std = std(ref_v)
-ref_v_a = ref_v[:-1]
-ref_v_b = ref_v[1:]
-ref_v_acf = pearsonr(ref_v_a, ref_v_b)
-ref_v_acf = ref_v_acf[0]
-
-# calculate theta rho
-acf = ref_rho_acf
-
-for i in range(len(xi)):
- if ia == 0 and xi[i+1] > -q:
- ia = i
- if ib == 0 and xi[i] > q:
- ib = i
- dnorm += 1 / sqrt(2*math.pi) * exp(-xi[i]*xi[i]/2)
-
-shape = ((NN+1)*(KK+1), 1)
-bb = matrix(zeros(shape))
-bb[-1:, 0] = dnorm
-
-Ta = [matrix(zeros((KK+1, KK+1))) for i in range(NN+1)]
-Tb = list(Ta)
-Tc = list(Ta)
-Td = [matrix(zeros((KK+1, 1))) for i in range(NN+1)]
-
-for i in range(NN+1):
- begin = (i)*(KK+1)
- end = begin + KK + 1
- Td[i] = bb[begin:end, 0]
-
-B1 = matrix(zeros((KK+1, KK+1)))
-B2 = matrix(B1)
-Id = matrix(B1)
-
-for i in range(KK+1):
- for j in range(KK+1):
- if i < ia or i >= ib:
- B1[i,j] = func_b(i,j)
- if i >= ia and i < ib:
- B2[i,j] = func_b(i,j)
- if i == j:
- Id[i,j] = 1
-
-Tb[0] = matrix(B2 - Id)
-Tc[0] = matrix(B2)
-for i in range(1, NN):
- Ta[i] = matrix(B1)
- Tb[i] = matrix(-Id)
- Tc[i] = matrix(B2)
-Ta[NN] = matrix(B1)
-Tb[NN] = matrix(B1 - Id)
-
-Tc[0] = linalg.solve(Tb[0], Tc[0])
-Td[0] = linalg.solve(Tb[0], Td[0])
-for i in range(1, NN+1):
- AA = Tb[i] - Ta[i] * Tc[i-1]
- Tc[i] = linalg.solve(AA, Tc[i])
- Td[i] = linalg.solve(AA, Td[i] - Ta[i] * Td[i-1])
-for i in range(NN-1, -1, -1):
- Td[i] = Td[i] - Tc[i] * Td[i+1]
-
-Tms = matrix(zeros(shape))
-for i in range(NN+1):
- begin = (i)*(KK+1)
- end = begin + KK + 1
- Tms[begin:end, 0] = Td[i]
-Tms = Tms / sum(d*Tms)
-
-shape = ((NN+1), 1)
-Tm = zeros(shape)
-for j in range(NN+1):
- begin = (j)*(KK+1)
- end = begin + KK + 1
- Tm[j, 0] = sum(d * Tms[begin:end, 0])
-
-Tps = Tm[0, 0]
-theta = 1
-while theta+1 < len(Tm) and Tps+Tm[theta, 0] < gamma:
- Tps = Tps + Tm[theta, 0]
- theta = theta + 1
-rho_theta = theta
-print('rho_theta = %.0f'%(rho_theta))
-
-# calculate theta v
-acf = ref_v_acf
-
-for i in range(len(xi)):
- if ia == 0 and xi[i+1] > -q:
- ia = i
- if ib == 0 and xi[i] > q:
- ib = i
- dnorm += 1 / sqrt(2*math.pi) * exp(-xi[i]*xi[i]/2)
-
-shape = ((NN+1)*(KK+1), 1)
-bb = matrix(zeros(shape))
-bb[-1:, 0] = dnorm
-
-Ta = [matrix(zeros((KK+1, KK+1))) for i in range(NN+1)]
-Tb = list(Ta)
-Tc = list(Ta)
-Td = [matrix(zeros((KK+1, 1))) for i in range(NN+1)]
-
-for i in range(NN+1):
- begin = (i)*(KK+1)
- end = begin + KK + 1
- Td[i] = bb[begin:end, 0]
-
-B1 = matrix(zeros((KK+1, KK+1)))
-B2 = matrix(B1)
-Id = matrix(B1)
-
-for i in range(KK+1):
- for j in range(KK+1):
- if i < ia or i >= ib:
- B1[i,j] = func_b(i,j)
- if i >= ia and i < ib:
- B2[i,j] = func_b(i,j)
- if i == j:
- Id[i,j] = 1
-
-Tb[0] = matrix(B2 - Id)
-Tc[0] = matrix(B2)
-for i in range(1, NN):
- Ta[i] = matrix(B1)
- Tb[i] = matrix(-Id)
- Tc[i] = matrix(B2)
-Ta[NN] = matrix(B1)
-Tb[NN] = matrix(B1 - Id)
-
-Tc[0] = linalg.solve(Tb[0], Tc[0])
-Td[0] = linalg.solve(Tb[0], Td[0])
-for i in range(1, NN+1):
- AA = Tb[i] - Ta[i] * Tc[i-1]
- Tc[i] = linalg.solve(AA, Tc[i])
- Td[i] = linalg.solve(AA, Td[i] - Ta[i] * Td[i-1])
-for i in range(NN-1, -1, -1):
- Td[i] = Td[i] - Tc[i] * Td[i+1]
-
-Tms = matrix(zeros(shape))
-for i in range(NN+1):
- begin = (i)*(KK+1)
- end = begin + KK + 1
- Tms[begin:end, 0] = Td[i]
-Tms = Tms / sum(d*Tms)
-
-shape = ((NN+1), 1)
-Tm = zeros(shape)
-for j in range(NN+1):
- begin = (j)*(KK+1)
- end = begin + KK + 1
- Tm[j, 0] = sum(d * Tms[begin:end, 0])
-
-Tps = Tm[0, 0]
-theta = 1
-while theta+1 < len(Tm) and Tps+Tm[theta, 0] < gamma:
- Tps = Tps + Tm[theta, 0]
- theta = theta + 1
-v_theta = theta
-print('v_theta = %.0f'%(v_theta))
-
-# calculate statistics rho
-statistics_rho = data
-statistics_rho = statistics_rho[:,1]
-file_rho_s = open('%s/cusum_rho_%s.txt'%(filepath, filename), 'w')
-file_rho_s.write('# frame s \n')
-rho_s_frame = 0
-rho_s = s_max
-file_rho_s.write('%.0f %.4f \n'%(rho_s_frame, rho_s))
-for i in statistics_rho:
- rho_s_frame = rho_s_frame + 1
- rho_s = min(max(0, rho_s + F(abs((i-ref_rho_mean)/ref_rho_std))), s_max)
- file_rho_s.write('%.0f %.4f \n'%(rho_s_frame, rho_s))
-file_rho_s.close()
-
-# calculate statistics v
-statistics_v = data
-statistics_v = statistics_v[:,2]
-file_v_s = open('%s/cusum_v_%s.txt'%(filepath, filename), 'w')
-file_v_s.write('# frame s \n')
-v_s_frame = 0
-v_s = s_max
-file_v_s.write('%.0f %.4f \n'%(v_s_frame, v_s))
-for i in statistics_v:
- v_s_frame = v_s_frame + 1
- v_s = min(max(0, v_s + F(abs((i-ref_v_mean)/ref_v_std))), s_max)
- file_v_s.write('%.0f %.4f \n'%(v_s_frame, v_s))
-file_v_s.close()
-
-# choose steady state rho
-statistics_rho = loadtxt('%s/cusum_rho_%s.txt'%(filepath, filename))
-steady_rho = statistics_rho[ statistics_rho[:,1] < rho_theta ]
-steady_rho_start_list = []
-steady_rho_end_list = []
-steady_rho_start = min(steady_rho[:,0]) - (s_max-rho_theta)
-steady_rho_start_list.extend([steady_rho_start])
-for i in arange(1, len(steady_rho), 1):
- if steady_rho[i,0] - steady_rho[i-1,0] != 1:
- steady_rho_end = steady_rho[i-1,0] - rho_theta
- steady_rho_end_list.extend([steady_rho_end])
- steady_rho_start = steady_rho[i,0] - (s_max-rho_theta)
- steady_rho_start_list.extend([steady_rho_start])
-if len(steady_rho_start_list) != len(steady_rho_end_list):
- steady_rho_end = max(steady_rho[:,0]) - rho_theta
- steady_rho_end_list.extend([steady_rho_end])
-steady_rho_list = []
-steady_rho_list.append(steady_rho_start_list)
-steady_rho_list.append(steady_rho_end_list)
-steady_rho_array = array(steady_rho_list)
-steady_rho_array = steady_rho_array.T
-savetxt('%s/SteadyState_rho_%s.txt'%(filepath, filename), steady_rho_array, fmt='%d\t%d', header='start end', newline='\r\n')
-print('steady state of rho is: \n', steady_rho_array)
-
-# choose steady state v
-statistics_v = loadtxt('%s/cusum_v_%s.txt'%(filepath, filename))
-steady_v = statistics_v[ statistics_v[:,1] < v_theta ]
-steady_v_start_list = []
-steady_v_end_list = []
-steady_v_start = min(steady_v[:,0]) - (s_max-v_theta)
-steady_v_start_list.extend([steady_v_start])
-for i in arange(1, len(steady_v), 1):
- if steady_v[i,0] - steady_v[i-1,0] != 1:
- steady_v_end = steady_v[i-1,0] - v_theta
- steady_v_end_list.extend([steady_v_end])
- steady_v_start = steady_v[i,0] - (s_max-v_theta)
- steady_v_start_list.extend([steady_v_start])
-if len(steady_v_start_list) != len(steady_v_end_list):
- steady_v_end = max(steady_v[:,0]) - v_theta
- steady_v_end_list.extend([steady_v_end])
-steady_v_list =[]
-steady_v_list.append(steady_v_start_list)
-steady_v_list.append(steady_v_end_list)
-steady_v_array = array(steady_v_list)
-steady_v_array = steady_v_array.T
-savetxt('%s/SteadyState_v_%s.txt'%(filepath, filename), steady_v_array, fmt='%d\t%d', header='start end', newline='\r\n')
-print('steady state of v is: \n', steady_v_array)
-
-# calculate steady state
-steady_start_list = []
-steady_end_list = []
-for i in range(len(steady_v_array[:,0])):
- for j in range(len(steady_rho_array[:,0])):
- if steady_v_array[i,0] <= steady_rho_array[j,0]:
- if steady_v_array[i,1] > steady_rho_array[j,0] and steady_v_array[i,1] < steady_rho_array[j,1]:
- steady_start = steady_rho_array[j,0]
- steady_end = steady_v_array[j,1]
- if steady_v_array[i,1] >= steady_rho_array[j,1]:
- steady_start = steady_rho_array[j,0]
- steady_end = steady_rho_array[j,1]
- if steady_v_array[i,0] > steady_rho_array[j,0] and steady_v_array[i,0] <= steady_rho_array[j,1]:
- if steady_v_array[i,1] <= steady_rho_array[j,1]:
- steady_start = steady_v_array[j,0]
- steady_end = steady_v_array[j,1]
- if steady_v_array[i,1] > steady_rho_array[j,1]:
- steady_start = steady_v_array[j,0]
- steady_end = steady_rho_array[j,1]
- steady_start_list.extend([steady_start])
- steady_end_list.extend([steady_end])
-steady_list =[]
-steady_list.append(steady_start_list)
-steady_list.append(steady_end_list)
-steady_array = array(steady_list)
-steady_array = steady_array.T
-savetxt('%s/SteadyState_%s.txt'%(filepath, filename), steady_array, fmt='%d\t%d', header='start end', newline='\r\n')
-print('final steady state is: \n', steady_array)
-
-if plotfigs == 'no':
- print('No figures are plotted!')
-else:
- print('Plotting figures...')
-
- # plot cusum rho
+def set_ref_series(data, column, ref_start, ref_end):
+ # calculate ref_mean, ref_std, ref_acf for data[:, column]
+ """
+ data:
+ column:
+ start:
+ end:
+ return: acf, mean, std
+ :param data:
+ :param column:
+ :param ref_start:
+ :param ref_end:
+ :return:
+ """
+ ref_series = data
+ ref_series = ref_series[ref_series[:, 0] >= ref_start - minframe]
+ ref_series = ref_series[ref_series[:, 0] <= ref_end - minframe]
+ ref_series = ref_series[:, column]
+ ref_series_mean = mean(ref_series)
+ ref_series_std = std(ref_series)
+ ref_series_a = ref_series[:-1]
+ ref_series_b = ref_series[1:]
+ ref_series_acf = pearsonr(ref_series_a, ref_series_b)
+ ref_series_acf = ref_series_acf[0]
+ return ref_series_acf, ref_series_mean, ref_series_std
+
+def get_theta(ia, ib, dnorm, acf):
+ """
+
+ :param ia:
+ :param ib:
+ :param dnorm:
+ :param acf:
+ :return:
+ """
+ shape = ((NN + 1) * (KK + 1), 1)
+ bb = matrix(zeros(shape))
+ bb[-1:, 0] = dnorm
+ Ta = [matrix(zeros((KK + 1, KK + 1))) for i in range(NN + 1)]
+ Tb = list(Ta)
+ Tc = list(Ta)
+ Td = [matrix(zeros((KK + 1, 1))) for i in range(NN + 1)]
+ for i in range(NN + 1):
+ begin = (i) * (KK + 1)
+ end = begin + KK + 1
+ Td[i] = bb[begin:end, 0]
+ B1 = matrix(zeros((KK + 1, KK + 1)))
+ B2 = matrix(B1)
+ Id = matrix(B1)
+ for i in range(KK + 1):
+ for j in range(KK + 1):
+ if i < ia or i >= ib:
+ B1[i, j] = func_b(i, j, acf)
+ if i >= ia and i < ib:
+ B2[i, j] = func_b(i, j, acf)
+ if i == j:
+ Id[i, j] = 1
+ Tb[0] = matrix(B2 - Id)
+ Tc[0] = matrix(B2)
+ for i in range(1, NN):
+ Ta[i] = matrix(B1)
+ Tb[i] = matrix(-Id)
+ Tc[i] = matrix(B2)
+ Ta[NN] = matrix(B1)
+ Tb[NN] = matrix(B1 - Id)
+ Tc[0] = linalg.solve(Tb[0], Tc[0])
+ Td[0] = linalg.solve(Tb[0], Td[0])
+ for i in range(1, NN + 1):
+ AA = Tb[i] - Ta[i] * Tc[i - 1]
+ Tc[i] = linalg.solve(AA, Tc[i])
+ Td[i] = linalg.solve(AA, Td[i] - Ta[i] * Td[i - 1])
+ for i in range(NN - 1, -1, -1):
+ Td[i] = Td[i] - Tc[i] * Td[i + 1]
+ Tms = matrix(zeros(shape))
+ for i in range(NN + 1):
+ begin = (i) * (KK + 1)
+ end = begin + KK + 1
+ Tms[begin:end, 0] = Td[i]
+ Tms = Tms / sum(d * Tms)
+ shape = ((NN + 1), 1)
+ Tm = zeros(shape)
+ for j in range(NN + 1):
+ begin = (j) * (KK + 1)
+ end = begin + KK + 1
+ Tm[j, 0] = sum(d * Tms[begin:end, 0])
+ Tps = Tm[0, 0]
+ theta = 1
+ while theta + 1 < len(Tm) and Tps + Tm[theta, 0] < gamma:
+ Tps = Tps + Tm[theta, 0]
+ theta += 1
+ rho_theta = theta
+ return rho_theta
+
+def init_parameters():
+ """
+ init ia, ib and dnorm
+ :return:
+ """
+ ia = ib = 0
+ dnorm = 0
+ for i in range(len(xi)):
+ if ia == 0 and xi[i + 1] > -q:
+ ia = i
+ if ib == 0 and xi[i] > q:
+ ib = i
+ dnorm += 1 / sqrt(2 * math.pi) * exp(-xi[i] * xi[i] / 2)
+
+ return ia, ib, dnorm
+
+def calculate_statistics(data, column, ref_mean, ref_std):
+ """
+ writes result in file
+ :param data:
+ :param column:
+ :param ref_mean:
+ :param ref_std:
+ """
+ statistics_series = data[:, column]
+ file_s = open('%s/cusum_%d_%s.txt' % (filepath, column, filename), 'w')
+ file_s.write('# frame s \n')
+ s_frame = minframe
+ s = s_max
+ file_s.write('%.0f %.4f \n' % (s_frame, s))
+ for i in statistics_series:
+ s_frame += 1
+ s = min(max(0, s + F(abs((i - ref_mean) / ref_std))), s_max)
+ file_s.write('%.0f %.4f \n' % (s_frame, s))
+ file_s.close()
+
+def choose_steady_state(column, theta):
+ """
+
+ :param column:
+ :param theta:
+ :return:
+ """
+ statistics = loadtxt('%s/cusum_%d_%s.txt' % (filepath, column, filename))
+ ss = open('%s/SteadyState_%d_%s.txt' % (filepath, column, filename), 'w')
+ ss.write('# start end ratio mean std \n')
+ steady = statistics[statistics[:, 1] < theta]
+ steady_start = min(steady[:, 0]) - (s_max - theta)
+ for i in arange(1, len(steady), 1):
+ if steady[i, 0] - steady[i - 1, 0] != 1:
+ steady_end = steady[i - 1, 0] - theta
+ if steady_start < steady_end:
+ series_data = data
+ series_data = series_data[series_data[:, 0] >= steady_start]
+ series_data = series_data[series_data[:, 0] <= steady_end]
+ data_ratio = len(series_data[:, 0]) / len(data[:, 0]) * 100
+ data_mean = mean(series_data[:, column])
+ data_std = std(series_data[:, column])
+ ss.write('%.0f %.0f %.2f %.4f %.4f \n' % (
+ steady_start, steady_end, data_ratio, data_mean, data_std))
+ steady_start = steady[i, 0] - (s_max - theta)
+ steady_end = max(steady[:, 0]) - theta
+ if steady_start < steady_end:
+ series_data = data
+ series_data = series_data[series_data[:, 0] >= steady_start]
+ series_data = series_data[series_data[:, 0] <= steady_end]
+ data_ratio = len(series_data[:, 0]) / len(data[:, 0]) * 100
+ data_mean = mean(series_data[:, column])
+ data_std = std(series_data[:, column])
+ ss.write('%.0f %.0f %.2f %.4f %.4f \n' % (
+ steady_start, steady_end, data_ratio, data_mean, data_std))
+ ss.close()
+ info_serie = loadtxt('%s/SteadyState_%d_%s.txt' % (filepath, column, filename))
+ if info_serie.shape == (5,):
+ temp = [info_serie]
+ info = array(temp)
+
+ return info, statistics
+
+def plot_series(statistics, theta, column):
fig = plt.figure(figsize=(11, 10), dpi=100)
- limit = (int((statistics_rho[-1,0]/frame)/10)+1)*10
- plt.plot(statistics_rho[:,0]/frame, statistics_rho[:,1], 'b--', lw=2, label=r'S$_{k}$')
- plt.plot([0,limit], [rho_theta, rho_theta], 'r-', lw=2, label=r'$\theta$')
- plt.xlabel('t [s]', fontsize=25)
- plt.ylabel('Statistics rho', fontsize=25)
+ limit = (int((statistics[-1, 0] / frame) / 10) + 1) * 10
+ plt.plot(statistics[:, 0] / frame, statistics[:, 1], 'b--', lw=2, label=r'S$_{k}$')
+ plt.plot([0, limit], [theta, theta], 'r-', lw=2, label=r'$\theta$')
+ plt.xlabel(xlabel, fontsize=25)
+ plt.ylabel('Statistics %d'%column, fontsize=25)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim(0, limit)
plt.ylim(-10, 120)
plt.legend(numpoints=2, ncol=1, loc=1, fontsize=20)
- plt.savefig('%s/cusum_rho_%s.png'%(filepath, filename))
-
- # plot cusum v
+ plt.savefig('%s/cusum_%d_%s.png' % (filepath, column, filename))
+ plt.close()
+
+
+def plot_steady_state(statistics, data, ref_start, ref_end, info, column):
fig = plt.figure(figsize=(11, 10), dpi=100)
- limit = (int((statistics_v[-1,0]/frame)/10)+1)*10
- plt.plot(statistics_v[:,0]/frame, statistics_v[:,1], 'b--', lw=2, label=r'S$_{k}$')
- plt.plot([0,limit], [v_theta, v_theta], 'r-', lw=2, label=r'$\theta$')
- plt.xlabel('t [s]', fontsize=25)
- plt.ylabel('Statistics v', fontsize=25)
- plt.xticks(fontsize=20)
- plt.yticks(fontsize=20)
- plt.xlim(0, limit)
- plt.ylim(-10, 120)
- plt.legend(numpoints=2, ncol=1, loc=1, fontsize=20)
- plt.savefig('%s/cusum_v_%s.png'%(filepath, filename))
-
- # plot steady rho
- fig = plt.figure(figsize=(11,10), dpi=100)
- ax = fig.add_subplot(111)
- limit = (int((statistics_rho[-1,0]/frame)/10)+1)*10
- plt.plot((data[:,0]-data[0,0])/frame, data[:,1], 'b-', lw=2)
- plt.plot([ref_rho_start/frame,ref_rho_start/frame], [0,50], 'g--', lw=2, label='reference process')
- plt.plot([ref_rho_end/frame,ref_rho_end/frame], [0,50], 'g--', lw=2)
- for i in range(len(steady_rho_array[:,0])):
- ax.add_patch(mpatches.Polygon([[steady_rho_array[i,0]/frame,0],[steady_rho_array[i,1]/frame,0],[steady_rho_array[i,1]/frame,50],[steady_rho_array[i,0]/frame,50]], closed=True, fill=False, color='r', hatch='/', label='steady state (rho)'))
- for i in range(len(steady_array[:,0])):
- ax.add_patch(mpatches.Polygon([[steady_array[i,0]/frame,0],[steady_array[i,1]/frame,0],[steady_array[i,1]/frame,50],[steady_array[i,0]/frame,50]], closed=True, fill=True, color='y', alpha=0.2, label='steady state (final)'))
- plt.xlabel('t [s]', fontsize=25)
- plt.ylabel(r'$\rho$ [m$^{-2}$]', fontsize=25)
- plt.xticks(fontsize=20)
- plt.yticks(fontsize=20)
- plt.xlim(0, limit)
- plt.ylim(0, int(max(data[:,1]))+2)
- plt.legend(numpoints=1, ncol=1, loc=1, fontsize=20)
- plt.savefig('%s/SteadyState_rho_%s.png'%(filepath, filename))
-
- # plot steady v
- fig = plt.figure(figsize=(11,10), dpi=100)
ax = fig.add_subplot(111)
- limit = (int((statistics_v[-1,0]/frame)/10)+1)*10
- plt.plot((data[:,0]-data[0,0])/frame, data[:,2], 'b-', lw=2)
- plt.plot([ref_v_start/frame,ref_v_start/frame], [0,50], 'g--', lw=2, label='reference process')
- plt.plot([ref_v_end/frame,ref_v_end/frame], [0,50], 'g--', lw=2)
- for i in range(len(steady_v_array[:,0])):
- ax.add_patch(mpatches.Polygon([[steady_v_array[i,0]/frame,0],[steady_v_array[i,1]/frame,0],[steady_v_array[i,1]/frame,50],[steady_v_array[i,0]/frame,50]], closed=True, fill=False, color='r', hatch='/', label='steady state (v)'))
- for i in range(len(steady_array[:,0])):
- ax.add_patch(mpatches.Polygon([[steady_array[i,0]/frame,0],[steady_array[i,1]/frame,0],[steady_array[i,1]/frame,50],[steady_array[i,0]/frame,50]], closed=True, fill=True, color='y', alpha=0.2, label='steady state (final)'))
- plt.xlabel('t [s]', fontsize=25)
- plt.ylabel(r'$v$ [m$\cdot$s$^{-1}$]', fontsize=25)
+ if xlabel == "frame":
+ fps = 1.0 # hack
+ else:
+ fps = frame
+
+ limit = (int((statistics[-1, 0] / fps) / 10) + 1) * 10
+ plt.plot((data[:, 0] + minframe) / fps, data[:, 1], 'b-', lw=2)
+ plt.plot([ref_start / fps, ref_start / fps], [0, 50], 'g--', lw=2, label='reference')
+ plt.plot([ref_end / fps, ref_end / fps], [0, 50], 'g--', lw=2)
+ for i in range(len(info[:, 0])):
+ ax.add_patch(mpatches.Polygon([[info[i, 0] / fps, 0],
+ [info[i, 1] / fps, 0],
+ [info[i, 1] / fps, 50],
+ [info[i, 0] / fps, 50]],
+ closed=True, fill=False,
+ color='r', hatch='/', label='steady (%d)'%column))
+
+ for i in range(len(info[:, 0])):
+ ax.add_patch(mpatches.Polygon([[info[i, 0] / fps, 0],
+ [info[i, 1] / fps, 0],
+ [info[i, 1] / fps, 50],
+ [info[i, 0] / fps, 50]],
+ closed=True, fill=True,
+ color='y', alpha=0.2, label='steady (final)'))
+
+ plt.xlabel(xlabel, fontsize=25)
+ plt.ylabel(r'$%s$'%ylabel[column-1], fontsize=25)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim(0, limit)
- plt.ylim(0, int(max(data[:,2]))+1)
+ plt.ylim(0, int(max(data[:, 1])) + 2)
plt.legend(numpoints=1, ncol=1, loc=1, fontsize=20)
- plt.savefig('%s/SteadyState_v_%s.png'%(filepath, filename))
+ plt.savefig('%s/SteadyState_%d_%s.png' % (filepath, column, filename))
+ plt.close()
+
+if __name__ == '__main__':
+ rho_max = 8.0
+ args = getParserArgs()
+ input_file = args.input_file
+ ref_start = args.reference_start
+ ref_end = args.reference_end
+ plotfigs = args.plotfigs
+ frame = args.fps
+ columns = args.columns
+ xlabel = args.xlabel
+ ylabel = args.ylabel
+ # sanity check
+
+ if not args.automatic: # in case references are manually given, lengths should be correct
+ assert (len(columns)-1) == len(ref_start) == len(ref_end) == len(ylabel),\
+ "mismatch lengths.\n\t columns: %s (first is frame)\n\t ref_start: %s\n\t ref_end: %s\n ylabel: %s"%\
+ (", ".join(map(str, columns)), ", ".join(map(str, ref_start)), ", ".join(map(str, ref_end)), ", ".join(map(str, ylabel)))
+
+ # read input data
+ try:
+ data = loadtxt('%s' % (input_file), usecols=columns)# [0, 6, 8]
+ except IOError:
+ exit("Can not open file <%s>" % input_file)
+
+ data = data[data[:, 1] != 0] # todo: why?
+ minframe = data[0, 0]
+ data[:, 0] = data[:, 0] - minframe
+
+ # get filepath and filename
+ filename = os.path.basename(input_file).split(".")[0]
+ filepath = os.path.dirname(input_file)
+ print('file path = %s' % filepath)
+ print('file name = %s' % filename)
+
+ ia, ib, dnorm = init_parameters()
+ starts = [] # collect start of steady state for all series
+ ends = [] # collect end of steady state for all series
+ for i in range(len(columns)-1):
+ ref_acf, ref_mean, ref_std = set_ref_series(data, i+1, ref_start[i], ref_end[i])
+ # calculate theta rho
+ theta = get_theta(ia, ib, dnorm, ref_acf)
+ print('theta[%d] = %.0f' % (i, theta))
+ # calculate statistics
+ calculate_statistics(data, i+1, ref_mean, ref_std)
+ # choose steady state
+ info, statistics = choose_steady_state(i+1, theta)
+ print('+------------------------------------------------------------------------------+')
+ for j in range(info.shape[0]):
+ print('steady state of series %d (%d): from %d (%.1f s) to %d (%.1f s) [ratio=%.2f, mean=%.2f, std=%.2f]' % (
+ i, j,
+ info[j][0], info[j][0] / frame, info[j][1], info[j][1] / frame,
+ info[j][2], info[j][3], info[j][4]))
+
+ ends.append(info[j][1])
+ starts.append(info[j][0])
+
+ if plotfigs == 'yes':
+ print('Plotting figures...')
+ # plot cusum
+ plot_series(statistics, theta, i+1)
+ # plot steady
+ plot_steady_state(statistics, data, ref_start[i], ref_end[i], info, i+1)
+
+ print('+------------------------------------------------------------------------------+')
+ os.remove('%s/cusum_%d_%s.txt' % (filepath, i+1, filename))
+
+ # choose steady state
+ ss = open('%s/SteadyState_%s.txt' % (filepath, filename), 'w')
+ ss.write('# start end ratio \n')
+ print "start frames: ", starts
+ print "end frames: ", ends
+ mix_start = max(starts)
+ mix_end = min(ends)
+ if mix_start < mix_end:
+ ss_data_ratio = (mix_end - mix_start) / len(data[:, 0]) * 100
+ ss.write('%.0f %.0f %.2f \n' % (mix_start, mix_end, ss_data_ratio))
+
+ ss.close()
+ info = loadtxt('%s/SteadyState_%s.txt' % (filepath, filename))
+ print('final steady state is from %d (%.1f s) to %d (%.1f s) [ratio=%.2f]' %
+ (mix_start, mix_start / frame, mix_end, mix_end / frame, ss_data_ratio))
+
+
+ print('Steady state detected successfully!')
+
+
+
+
+ # choose steady state v
+ # info_v, statistics_v = choose_steady_state(2, v_theta)
+
+ # print('+--------------------------------------------------------------------------------------------+')
+ # for i in range(info_v.shape[0]):
+ # print('steady state of v (%d): from %d (%.1f s) to %d (%.1f s) [ratio=%.2f, mean=%.2f, std=%.2f]' % (
+ # i,
+ # info_v[i][0], info_v[i][0] / frame, info_v[i][1], info_v[i][1] / frame,
+ # info_v[i][2], info_v[i][3], info_v[i][4]))
+ # print('+--------------------------------------------------------------------------------------------+')
+
+ # calculate steady state
+
+
+ # for i in range(len(info_rho[:, 0])):
+ # for j in range(len(info_v[:, 0])):
+ # mix_start = max(info_rho[i, 0], info_v[j, 0])
+ # mix_end = min(info_rho[i, 1], info_v[j, 1])
+ # if mix_start < mix_end:
+ # ss_data_ratio = (mix_end - mix_start) / len(data[:, 0]) * 100
+ # ss.write('%.0f %.0f %.2f \n' % (mix_start, mix_end, ss_data_ratio))
+ # ss.close()
+ # info = loadtxt('%s/SteadyState_%s.txt' % (filepath, filename))
+ # if info.shape == (3,):
+ # temp = [info]
+ # info = array(temp)
+
+ # print('final steady state is from %d (%.1f s) to %d (%.1f s) [ratio=%.2f]' %
+ # (info[0][0], info[0][0] / frame, info[0][1], info[0][1] / frame, info[0][2]))
+ # print('+--------------------------------------------------------------------------------------------+')
-print('Steady state detected successfully!')
\ No newline at end of file
diff --git a/scripts/_Plot_FD.py b/scripts/_Plot_FD.py
index f786d81a1023155de9327533da0b2077aeee05da..c0269c1bdafe74b74eb2099fcf5e86b9510cf9ce 100644
--- a/scripts/_Plot_FD.py
+++ b/scripts/_Plot_FD.py
@@ -29,7 +29,7 @@ def plotRhoJs(pathfile,figname,namefile,steadyfile):
steady_rho_v=zeros(len(d_rho_v[0]))
for i in range(len(d_steady)):
rho_v=d_rho_v[d_rho_v[:,0]>d_steady[i,0]]
- rho_v=rho_v[rho_v[:,0]>d_steady[i,1]]
+ rho_v=rho_v[rho_v[:,0]<d_steady[i,1]]
steady_rho_v=vstack((steady_rho_v,rho_v))
steady_rho_v = delete(steady_rho_v, (0), axis=0)
plt.plot(steady_rho_v[:,1],steady_rho_v[:,1]*steady_rho_v[:,2], 'o', lw=3, label=f_rho_v)
@@ -53,7 +53,7 @@ def plotRhoV(pathfile,figname,namefile,steadyfile):
steady_rho_v=zeros(len(d_rho_v[0]))
for i in range(len(d_steady)):
rho_v=d_rho_v[d_rho_v[:,0]>d_steady[i,0]]
- rho_v=rho_v[rho_v[:,0]>d_steady[i,1]]
+ rho_v=rho_v[rho_v[:,0]<d_steady[i,1]]
steady_rho_v=vstack((steady_rho_v,rho_v))
steady_rho_v = delete(steady_rho_v, (0), axis=0)
plt.plot(steady_rho_v[:,1],steady_rho_v[:,2], 'o', lw=3, label=f_rho_v)
diff --git a/scripts/_Plot_cell_rho.py b/scripts/_Plot_cell_rho.py
index f55b74441d49be7d1da1b6074155d0c8ebec90da..437dbecbfb26fd28a308cae315b8ae45b736410f 100644
--- a/scripts/_Plot_cell_rho.py
+++ b/scripts/_Plot_cell_rho.py
@@ -8,19 +8,135 @@ import pylab
from mpl_toolkits.axes_grid1 import make_axes_locatable
import argparse
import sys
-
+import logging
+from xml.dom import minidom
+import os
+try:
+ import xml.etree.cElementTree as ET
+except ImportError:
+ import xml.etree.ElementTree as ET
+
def getParserArgs():
parser = argparse.ArgumentParser(description='Combine French data to one file')
parser.add_argument("-f", "--filepath", default="./", help='give the path of source file')
parser.add_argument("-n", "--namefile", help='give the name of the source file')
- parser.add_argument("-x1", "--geominx", type=float, help='give the minmum x of the geometry')
- parser.add_argument("-x2", "--geomaxx", type=float, help='give the maxmum x of the geometry')
- parser.add_argument("-y1", "--geominy", type=float, help='give the minmum y of the geometry')
- parser.add_argument("-y2", "--geomaxy", type=float, help='give the maxmum y of the geometry')
+ parser.add_argument("-g", "--geoname", help='give the name of the geometry file')
args = parser.parse_args()
return args
+def get_polygon(poly):
+ X = []
+ Y = []
+ # for poly in node.getchildren():
+ # # print "polygon tag: ", poly.tag
+ # # if poly.tag == "obstacle":
+ # # # print poly.getchildren()
+ # # # for pobst in poly.getchildren():
+ # # # # print pobst.tag
+ # # # for q in pobst.getchildren(): # vertex
+ # # # X.append( q.attrib['px'] )
+ # # # Y.append( q.attrib['py'] )
+ # # pass
+ # # else:
+ for p in poly.getchildren(): # vertex
+ X.append(float(p.attrib['px']))
+ Y.append(float(p.attrib['py']))
+
+ return X, Y
+
+def get_geometry_boundary(geometry):
+ tree = ET.parse(geometry)
+ root = tree.getroot()
+ geominX = []
+ geomaxX = []
+ geominY = []
+ geomaxY = []
+ for node in root.iter():
+ tag = node.tag
+ # print "subroom tag", tag
+ if tag == "polygon":
+ X, Y = get_polygon(node)
+ geominX.append(min(X))
+ geomaxX.append(max(X))
+ geominY.append(min(Y))
+ geomaxY.append(max(Y))
+ elif tag == "obstacle":
+ # print "obstacle tag",tag
+ for n in node.getchildren():
+ X, Y = get_polygon(n)
+ geominX.append(min(X))
+ geomaxX.append(max(X))
+ geominY.append(min(Y))
+ geomaxY.append(max(Y))
+ geominX = min(geominX)
+ geomaxX = max(geomaxX)
+ geominY = min(geominY)
+ geomaxY = max(geomaxY)
+ return geominX, geomaxX, geominY, geomaxY
+
+def plot_geometry(geometry):
+ tree = ET.parse(geometry)
+ root = tree.getroot()
+ for node in root.iter():
+ tag = node.tag
+ # print "subroom tag", tag
+ if tag == "polygon":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "k", lw=4)
+ elif tag == "obstacle":
+ # print "obstacle tag",tag
+ for n in node.getchildren():
+ # print "N: ", n
+ X, Y = get_polygon(n)
+ # print "obstacle", X, Y
+ plt.plot(X, Y, "g", lw=4)
+ elif tag == "crossing":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "--b", lw=0.9, alpha=0.2)
+ elif tag == "HLine":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "--b", lw=0.9, alpha=0.2)
+ elif tag == "transition":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "--r", lw=1.6, alpha=0.2)
+
+
+def parse_xml_traj_file(filename):
+ """
+ parse trajectories in Travisto-format and output results
+ in the following format: id frame x y
+ (no sorting of the data is performed)
+ returns:
+ fps: frames per second
+ N: number of pedestrians
+ data: trajectories (numpy.array) [id fr x y]
+ """
+ logging.info("parsing <%s>"%filename)
+ try:
+ xmldoc = minidom.parse(filename)
+ except:
+ logging.critical('could not parse file. exit')
+ exit(FAILURE)
+ N = int(xmldoc.getElementsByTagName('agents')[0].childNodes[0].data)
+ fps = xmldoc.getElementsByTagName('frameRate')[0].childNodes[0].data #type unicode
+ fps = float(fps)
+ fps = int(fps)
+ #print ("fps=", fps)
+ #fps = int(xmldoc.getElementsByTagName('frameRate')[0].childNodes[0].data)
+ logging.info("Npeds = %d, fps = %d"%(N, fps))
+ frames = xmldoc.childNodes[0].getElementsByTagName('frame')
+ data = []
+ for frame in frames:
+ frame_number = int(frame.attributes["ID"].value)
+ for agent in frame.getElementsByTagName("agent"):
+ agent_id = int(agent.attributes["ID"].value)
+ x = float(agent.attributes["x"].value)
+ y = float(agent.attributes["y"].value)
+ data += [agent_id, frame_number, x, y]
+ data = array(data).reshape((-1, 4))
+ return fps, N, data
+
if __name__ == '__main__':
rho_max = 8.0
@@ -28,10 +144,14 @@ if __name__ == '__main__':
filepath = args.filepath
sys.path.append(filepath)
namefile = args.namefile
- geominX = args.geominx
- geomaxX = args.geomaxx
- geominY = args.geominy
- geomaxY = args.geomaxy
+ geoFile=args.geoname
+ geoLocation = filepath.split("Output")[0]
+ trajName = namefile.split(".")[0]
+ trajType = namefile.split(".")[1].split("_")[0]
+ trajFile = geoLocation+trajName+"."+trajType
+ frameNr=int(namefile.split("_")[-1])
+ geominX, geomaxX, geominY, geomaxY = get_geometry_boundary(geoFile)
+
fig = plt.figure(figsize=(16*(geomaxX-geominX)/(geomaxY-geominY)+2, 16), dpi=100)
ax1 = fig.add_subplot(111,aspect='equal')
plt.rc("font", size=30)
@@ -42,12 +162,13 @@ if __name__ == '__main__':
for tick in ax1.get_xticklines() + ax1.get_yticklines():
tick.set_markeredgewidth(2)
tick.set_markersize(6)
- ax1.set_aspect("equal")
+ ax1.set_aspect("equal")
+ plot_geometry(geoFile)
+
density=array([])
polys = open("%s/polygon%s.dat"%(filepath,namefile)).readlines()
for poly in polys:
exec("p = %s"%poly)
- #p=tuple(tuple(i/100.0 for i in inner) for inner in p)
xx=1.0/Polygon(p).area()
if xx>rho_max:
xx=rho_max
@@ -62,19 +183,22 @@ if __name__ == '__main__':
sm.set_clim(vmin=0,vmax=5)
for poly in polys:
exec("p = %s"%poly)
- #p=tuple(tuple(i/100.0 for i in inner) for inner in p)
xx=1.0/Polygon(p).area()
- #print xx
if xx>rho_max:
xx=rho_max
ax1.add_patch(pgon(p,facecolor=sm.to_rgba(xx), edgecolor='white',linewidth=2))
- points = loadtxt("%s/points%s.dat"%(filepath,namefile))
- ax1.plot(points[:,0],points[:,1],"bo",markersize = 20,markeredgewidth=2)
- ax1.set_xlim(geominX,geomaxX)
- ax1.set_ylim(geominY,geomaxY)
+
+ if(trajType=="xml"):
+ fps, N, Trajdata = parse_xml_traj_file(trajFile)
+ elif(trajType=="txt"):
+ Trajdata = loadtxt(trajFile)
+ Trajdata = Trajdata[ Trajdata[:, 1] == frameNr ]
+ ax1.plot(Trajdata[:,2], Trajdata[:,3], "bo", markersize = 20, markeredgewidth=2)
+
+ ax1.set_xlim(geominX-0.2,geomaxX+0.2)
+ ax1.set_ylim(geominY-0.2,geomaxY+0.2)
plt.xlabel("x [m]")
plt.ylabel("y [m]")
- #print density
plt.title("%s"%namefile)
divider = make_axes_locatable(ax1)
cax = divider.append_axes("right", size="2.5%", pad=0.2)
diff --git a/scripts/_Plot_cell_v.py b/scripts/_Plot_cell_v.py
index a9851d75c66e52ad67e703b010be2b0348144ddc..74c1ffce3fd10ae2c71c300427e4265d6b9dca0b 100644
--- a/scripts/_Plot_cell_v.py
+++ b/scripts/_Plot_cell_v.py
@@ -8,16 +8,14 @@ from mpl_toolkits.axes_grid1 import make_axes_locatable
import matplotlib
import argparse
import sys
+from _Plot_cell_rho import *
def getParserArgs():
parser = argparse.ArgumentParser(description='Combine French data to one file')
parser.add_argument("-f", "--filepath", default="./", help='give the path of source file')
parser.add_argument("-n", "--namefile", help='give the name of the source file')
- parser.add_argument("-x1", "--geominx", type=float, help='give the minmum x of the geometry')
- parser.add_argument("-x2", "--geomaxx", type=float, help='give the maxmum x of the geometry')
- parser.add_argument("-y1", "--geominy", type=float, help='give the minmum y of the geometry')
- parser.add_argument("-y2", "--geomaxy", type=float, help='give the maxmum y of the geometry')
+ parser.add_argument("-g", "--geoname", help='give the name of the geometry file')
args = parser.parse_args()
return args
@@ -27,10 +25,14 @@ if __name__ == '__main__':
filepath = args.filepath
sys.path.append(filepath)
namefile = args.namefile
- geominX = args.geominx
- geomaxX = args.geomaxx
- geominY = args.geominy
- geomaxY = args.geomaxy
+ geoFile=args.geoname
+ geoLocation = filepath.split("Output")[0]
+ trajName = namefile.split(".")[0]
+ trajType = namefile.split(".")[1].split("_")[0]
+ trajFile = geoLocation+trajName+"."+trajType
+ frameNr=int(namefile.split("_")[-1])
+ geominX, geomaxX, geominY, geomaxY = get_geometry_boundary(geoFile)
+
fig = plt.figure(figsize=(16*(geomaxX-geominX)/(geomaxY-geominY)+2, 16), dpi=100)
ax1 = fig.add_subplot(111,aspect='equal')
plt.rc("font", size=30)
@@ -40,7 +42,9 @@ if __name__ == '__main__':
for tick in ax1.get_xticklines() + ax1.get_yticklines():
tick.set_markeredgewidth(2)
tick.set_markersize(6)
- ax1.set_aspect("equal")
+ ax1.set_aspect("equal")
+ plot_geometry(geoFile)
+
velocity=array([])
polys = open("%s/polygon%s.dat"%(filepath,namefile)).readlines()
velocity=loadtxt("%s/speed%s.dat"%(filepath,namefile))
@@ -51,14 +55,19 @@ if __name__ == '__main__':
index=0
for poly in polys:
exec("p = %s"%poly)
- #p=tuple(tuple(i/100.0 for i in inner) for inner in p)
xx = velocity[index]
index += 1
ax1.add_patch(pgon(p,facecolor=sm.to_rgba(xx), edgecolor='white',linewidth=2))
- points = loadtxt("%s/points%s.dat"%(filepath,namefile))
- ax1.plot(points[:,0],points[:,1],"bo",markersize = 20,markeredgewidth=2)
- ax1.set_xlim(geominX,geomaxX)
- ax1.set_ylim(geominY,geomaxY)
+
+ if(trajType=="xml"):
+ fps, N, Trajdata = parse_xml_traj_file(trajFile)
+ elif(trajType=="txt"):
+ Trajdata = loadtxt(trajFile)
+ Trajdata = Trajdata[ Trajdata[:, 1] == frameNr ]
+ ax1.plot(Trajdata[:,2], Trajdata[:,3], "bo", markersize = 20, markeredgewidth=2)
+
+ ax1.set_xlim(geominX-0.2,geomaxX+0.2)
+ ax1.set_ylim(geominY-0.2,geomaxY+0.2)
plt.xlabel("x [m]")
plt.ylabel("y [m]")
plt.title("%s"%namefile)
diff --git a/scripts/_Plot_timeseries_rho_v.py b/scripts/_Plot_timeseries_rho_v.py
index cb3b68b96220dc0e3eebf34d385179512e26f10e..b296123ef94e04a544d8cdbb335c0bc5e5a7d047 100644
--- a/scripts/_Plot_timeseries_rho_v.py
+++ b/scripts/_Plot_timeseries_rho_v.py
@@ -86,5 +86,5 @@ if __name__ == '__main__':
plotRhoT(pathfile,figname_rho,fps,title,data_Classic,data_Voronoi)
plotVT(pathfile,figname_v,fps,title,data_Classic,data_Voronoi)
else:
- plotRhoT(pathfile,figname_rho,fps,title,data_Voronoi)
- plotVT(pathfile,figname_v,fps,title,data_Voronoi)
+ plotRhoT(pathfile,figname_rho,fps,title,data_Voronoi=data_Voronoi)
+ plotVT(pathfile,figname_v,fps,title,data_Voronoi=data_Voronoi)
diff --git a/scripts/plot_geom.py b/scripts/plot_geom.py
new file mode 100644
index 0000000000000000000000000000000000000000..75a1b4d45def19858c19d6810d2c53b3bfedbf6a
--- /dev/null
+++ b/scripts/plot_geom.py
@@ -0,0 +1,62 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from sys import argv
+
+try:
+ import xml.etree.cElementTree as ET
+except ImportError:
+ import xml.etree.ElementTree as ET
+
+
+def get_polygon(poly):
+ X = []
+ Y = []
+ # for poly in node.getchildren():
+ # # print "polygon tag: ", poly.tag
+ # # if poly.tag == "obstacle":
+ # # # print poly.getchildren()
+ # # # for pobst in poly.getchildren():
+ # # # # print pobst.tag
+ # # # for q in pobst.getchildren(): # vertex
+ # # # X.append( q.attrib['px'] )
+ # # # Y.append( q.attrib['py'] )
+ # # pass
+ # # else:
+ for p in poly.getchildren(): # vertex
+ X.append(p.attrib['px'])
+ Y.append(p.attrib['py'])
+
+ return X, Y
+
+def plot_geometry(filename):
+ tree = ET.parse(geometry)
+ root = tree.getroot()
+ for node in root.iter():
+ tag = node.tag
+ # print "subroom tag", tag
+ if tag == "polygon":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "k", lw=2)
+ elif tag == "obstacle":
+ # print "obstacle tag",tag
+ for n in node.getchildren():
+ # print "N: ", n
+ X, Y = get_polygon(n)
+ # print "obstacle", X, Y
+ plt.plot(X, Y, "g", lw=2)
+ elif tag == "crossing":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "--b", lw=0.9, alpha=0.2)
+ elif tag == "HLine":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "--b", lw=0.9, alpha=0.2)
+ elif tag == "transition":
+ X, Y = get_polygon(node)
+ plt.plot(X, Y, "--r", lw=1.6, alpha=0.2)
+
+
+
+geometry = argv[1]
+plot_geometry(geometry)
+
+plt.show()