Commit a0e7970d authored by Mohcine Chraibi's avatar Mohcine Chraibi

Hack for the polar radius in wallRepForce

- Assume a constant value (actually 0.4 m) but could be something
  between bmax and amin + v*atau
parent 4d81b558
......@@ -90,3 +90,18 @@ vgcore.*
/delme.py
/mencoderVideos.sh
/true_false.xml
/Max/
README.html
Utest/r2.txt
Utest/r3.txt
Utest/rimea_tests/test_13/gauss.py
Utest/test_10/cell.png
Utest/test_10/result.txt
Utest/test_10/times.txt
Utest/todo.org
html-dir/
inputfiles/Bottleneck/rimea_test_12_geo.xml
inputfiles/Karthik/
inputfiles/U8/
inputfiles/VerteilerEbene/
......@@ -395,10 +395,20 @@ if(CMAKE_COMPILER_IS_GNUCXX OR CMAKE_CXX_COMPILER_ID STREQUAL "Clang")
endif()
#--------------------
if(APPLE)
#------
SET( CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -stdlib=libc++ -std=c++11" )
SET( CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} -stdlib=libc++" )
SET( CMAKE_SHARED_LINKER_FLAGS "${CMAKE_SHARED_LINKER_FLAGS} -stdlib=libc++" )
SET( CMAKE_MODULE_LINKER_FLAGS "${CMAKE_MODULE_LINKER_FLAGS} -stdlib=libc++" )
#------
endif()
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)
if( COMPILER_SUPPORTS_CXX11)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
elseif(COMPILER_SUPPORTS_CXX0X)
......@@ -407,6 +417,7 @@ if(NOT CMAKE_GENERATOR MATCHES "Xcode|Visual Studio")
message(FATAL_ERROR "The compiler ${CMAKE_CXX_COMPILER} has no C++11 support. Please use a different C++ compiler.")
endif()
endif()
message(STATUS "MSVC: " ${MSVC})
message(STATUS "MSVC_IDE: " ${MSVC_IDE})
message(STATUS "MSVC60: " ${MSVC60})
......
......@@ -2,7 +2,7 @@
<!-- <JuPedSim project="JPS-Project" version="0.5" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:noNamespaceSchemaLocation="ini.xsd"> -->
<JuPedSim project="JPS-Project" version="0.5" xmlns:xsi="http://xsd.jupedsim.org/jps_ini_core.xsd" xsi:noNamespaceSchemaLocation="ini.xsd">
<seed>linspace(1,10000,3)</seed>
<seed>linspace(1,10000,1)</seed>
<!-- <seed>1</seed> -->
<num_threads>4</num_threads>
<max_sim_time unit="sec">100</max_sim_time>
......@@ -30,15 +30,15 @@
<routing>
<goals>
<goal id="0" final="true" caption="goal 0">
<polygon>
<vertex px="11.0" py="2.0" />
<vertex px="11.0" py="3.0" />
<vertex px="12.0" py="3.0" />
<vertex px="12.0" py="2.0" />
<vertex px="11.0" py="2.0" />
</polygon>
</goal>
<!-- <goal id="0" final="true" caption="goal 0"> -->
<!-- <polygon> -->
<!-- <vertex px="11.0" py="2.0" /> -->
<!-- <vertex px="11.0" py="3.0" /> -->
<!-- <vertex px="12.0" py="3.0" /> -->
<!-- <vertex px="12.0" py="2.0" /> -->
<!-- <vertex px="11.0" py="2.0" /> -->
<!-- </polygon> -->
<!-- </goal> -->
</goals>
</routing>
......@@ -58,7 +58,7 @@
<model_parameters>
<solver>euler</solver>
<stepsize>0.001</stepsize>
<exit_crossing_strategy>4</exit_crossing_strategy>
<exit_crossing_strategy>2</exit_crossing_strategy>
<linkedcells enabled="true" cell_size="2.2" />
<force_ped nu="0.2" dist_max="2" disteff_max="2" interpolation_width="0.1" />
<force_wall nu="0.2" dist_max="1" disteff_max="1" interpolation_width="0.1" />
......@@ -79,8 +79,11 @@
<stepsize>0.001</stepsize>
<exit_crossing_strategy>3</exit_crossing_strategy>
<linkedcells enabled="true" cell_size="3" />
<force_ped nu="0.5" b="0.06" c="3.45"/>
<force_wall nu="5" b="0.4" c="10.0"/>
<!-- <force_ped nu="0.5" b="0.06" c="3.45"/> -->
<!-- <force_wall nu="5" b="0.4" c="10.0"/> -->
<force_ped nu="0." b="0.034229" c="2.450932"/>
<force_wall nu="0." b="0.034229" c="2.450932"/>
</model_parameters>
<agent_parameters agent_parameter_id="0">
<v0 mu="1.34" sigma="0.0" />
......
......@@ -140,10 +140,8 @@ def flow(fps, N, data, x0):
if len(times) < 2:
logging.warning("Number of pedestrians passing the line is small. return 0")
return 0
logging.info("min(times)=%f max(times)=%f"%(min(times)/fps, max(times)/fps))
flow = fps * float(N-1) / (max(times) - min(times))
return flow
logging.info("min(times)=%f max(times)=%f"%(min(times), max(times)))
return fps * float(N-1) / (max(times) - min(times))
......
......@@ -188,24 +188,25 @@ Point Line::LotPoint(const Point &p) const {
return f;
}
/* Punkt auf der Linie mit kürzestem Abstand zu p
* In der Regel Lotfußpunkt, Ist der Lotfußpunkt nicht im Segment
* wird der entsprechende Eckpunkt der Line genommen
* */
Point Line::ShortestPoint(const Point &p) const {
if (_point1 == _point2)
return _point1;
const Point &t = _point1 - _point2;
double lambda = (p - _point2).ScalarProduct(t) / t.ScalarProduct(t);
if(lambda <0)
return _point2;
else if(lambda >1)
return _point1;
else
return _point2 + t * lambda;
}
Point Line::ShortestPoint(const Point &p) const {
if (_point1 == _point2)
return _point1;
const Point &t = _point1 - _point2;
double lambda = (p - _point2).ScalarProduct(t) / t.ScalarProduct(t);
if(lambda <0)
return _point2;
else if(lambda >1)
return _point1;
else
return _point2 + t * lambda;
}
/*
* Prüft, ob Punkt p im Liniensegment enthalten ist
......
......@@ -8,7 +8,7 @@
<seed>12542</seed>
<max_sim_time >900</max_sim_time>
<!-- geometry file -->
<geometry>2.0_bottleneck.xml</geometry>
<geometry>1.2_bottleneck.xml</geometry>
<!-- traectories file and format -->
<trajectories format="xml-plain" fps="8">
<file location="trajectorien_bottleneck.xml" />
......@@ -47,7 +47,10 @@
<!--persons information and distribution -->
<agents operational_model_id="2">
<agents_distribution>
<group group_id="2" agent_parameter_id="1" room_id="1" subroom_id="0" number="10" goal_id="0" router_id="1" />
<group group_id="2" agent_parameter_id="1" room_id="1" subroom_id="0" number="3
pwd
0" goal_id="0" router_id="1" />
<!-- <group group_id="0" room_id="0" subroom_id="2" number="4" goal_id="0" router_id="1" route_id="" /> -->
<!-- <group group_id="1" room_id="1" number="10" goal_id="0" router_id="1" route_id="" /> -->
</agents_distribution>
......@@ -82,8 +85,8 @@
<stepsize>0.001</stepsize>
<exit_crossing_strategy>4</exit_crossing_strategy>
<linkedcells enabled="true" cell_size="2.2" />
<force_ped nu="1" b="0.06" c="3.45"/>
<force_wall nu="5" b="0.4" c="10"/>
<force_ped nu="0.5" b="0.034229" c="2.450932"/>
<force_wall nu="1" b="0.034229" c="2.450932"/>
</model_parameters>
<agent_parameters agent_parameter_id="1">
<v0 mu="1.34" sigma="0.001" />
......
......@@ -193,6 +193,12 @@ void GompertzModel::ComputeNextTimeStep(double current, double deltaT, Building*
Point fd = ForceDriv(ped, room);
Point acc = (fd + repPed + repWall) / ped->GetMass();
if(ped->GetID()==-4)
{
printf("t=%f, Pos1 =[%f, %f]\n", current,ped->GetPos().GetX(), ped->GetPos().GetY());
printf("acc= %f %f, fd= %f, %f, repPed = %f %f, repWall= %f, %f\n", acc.GetX(), acc.GetY(), fd.GetX(), fd.GetY(), repPed.GetX(), repPed.GetY(), repWall.GetX(), repWall.GetY());
// if(current >16) getc(stdin);
}
result_acc.push_back(acc);
}
......@@ -247,23 +253,28 @@ Point GompertzModel::ForceDriv(Pedestrian* ped, Room* room) const
// check if the molified version works
if (dist > J_EPS_GOAL) {
e0 = ped->GetV0(target);
// printf("1 e0 %f %f, target %f %f\n", e0.GetX(), e0.GetY(), target.GetX(), target.GetY());
if(ped->GetID()==-4)
printf("1 e0 %f %f, target %f %f\n", e0.GetX(), e0.GetY(), target.GetX(), target.GetY());
} else {
ped->SetSmoothTurning();
e0 = ped->GetV0();
// printf("2 e0 %f %f\n", e0.GetX(), e0.GetY());
if(ped->GetID()==-4)
printf("2 e0 %f %f\n", e0.GetX(), e0.GetY());
}
F_driv = ((e0 * ped->GetV0Norm() - ped->GetV()) * ped->GetMass()) / ped->GetTau();
//double v = sqrt(ped->GetV().GetX()*ped->GetV().GetX() +ped->GetV().GetY()*ped->GetV().GetY());
#if DEBUG
double e0norm = sqrt(e0.GetX()*e0.GetX() +e0.GetY()*e0.GetY());
printf( "pos %f %f target %f %f\n", pos.GetX(), pos.GetY(), target.GetX(), target.GetY());
printf("mass=%f, v0norm=%f, e0Norm=%f, tau=%f\n", ped->GetMass(), ped->GetV0Norm(), e0norm,ped->GetTau());
printf("Fdriv= [%f, %f]\n", F_driv.GetX(), F_driv.GetY());
fprintf(stderr, "%d %f %f %f %f %f %f\n", ped->GetID(), ped->GetPos().GetX(), ped->GetPos().GetY(), ped->GetV().GetX(), ped->GetV().GetY(), target.GetX(), target.GetY());
#endif
// #if DEBUG
if (ped->GetID()==-4){
double e0norm = sqrt(e0.GetX()*e0.GetX() +e0.GetY()*e0.GetY());
printf( "pos %f %f target %f %f\n", pos.GetX(), pos.GetY(), target.GetX(), target.GetY());
printf("mass=%f, v0norm=%f, e0Norm=%f, tau=%f\n", ped->GetMass(), ped->GetV0Norm(), e0norm, ped->GetTau());
printf("Fdriv= [%f, %f]\n", F_driv.GetX(), F_driv.GetY());
fprintf(stdout, "%d %f %f %f %f %f %f\n", ped->GetID(), ped->GetPos().GetX(), ped->GetPos().GetY(), ped->GetV().GetX(), ped->GetV().GetY(), target.GetX(), target.GetY());
}
// #endif
return F_driv;
}
......@@ -309,7 +320,7 @@ Point GompertzModel::ForceRepPed(Pedestrian* ped1, Pedestrian* ped2) const
}
//------------------------- check if others are behind using v0 instead of v
double tmpv = ped1->GetV().ScalarProduct(ep12); // < v^0_i , e_ij >
double ped2IsBehindv = exp(-exp(-5*tmpv)); //step function: continuous version
double ped2IsBehindv = (tmpv<=0)?0:1; //exp(-exp(-5*tmpv)); //step function: continuous version
if (ped2IsBehindv < J_EPS) {
return F_rep; // ignore ped2
}
......@@ -327,9 +338,13 @@ Point GompertzModel::ForceRepPed(Pedestrian* ped1, Pedestrian* ped2) const
f = -ped1->GetMass() * _nuPed * ped1->GetV0Norm() * B_ij;
F_rep = ep12 * f;
// if(ped1->GetID() == 1) {
// printf("F=[%f, %f] v0=%f, nu=%f, B_ij=%f D=%f, r1=%f, r2=%f\n", F_rep.GetX(), F_rep.GetY(), ped1->GetV0Norm(), _nuPed, B_ij, Distance, r1, r2);
// }
if(ped1->GetID() ==-4) {
printf("\nNAN return ----> p1=%d p2=%d pos1=%f %f, pos2=%f %f\n", ped1->GetID(),
ped2->GetID(), ped1->GetPos().GetX(), ped1->GetPos().GetY(), ped2->GetPos().GetX(), ped2->GetPos().GetY());
printf("F=[%f, %f] v0=%f, nu=%f, B_ij=%f D=%f, r1=%f, r2=%f\n", F_rep.GetX(), F_rep.GetY(), ped1->GetV0Norm(), _nuPed, B_ij, Distance, r1, r2);
}
//check isNan
if (F_rep.GetX() != F_rep.GetX() || F_rep.GetY() != F_rep.GetY()) {
char tmp[CLENGTH];
......@@ -395,13 +410,14 @@ Point GompertzModel::ForceRepWall(Pedestrian* ped, const Line& w, const Point& c
#define DEBUG 0
Point F_wrep = Point(0.0, 0.0);
#if DEBUG
if(ped->GetID()==-4)
printf("=========\n\tEnter GompertzWall with PED=%d, wall=[%.2f, %.2f]--[%.2f, %.2f]\n", ped->GetID(), w.GetPoint1().GetX(), w.GetPoint1().GetY(), w.GetPoint2().GetX(), w.GetPoint2().GetY());
#endif
// getc(stdin);
// if direction of pedestrians does not intersect walls --> ignore
Point pt = w.ShortestPoint(ped->GetPos());
double wlen = w.LengthSquare();
// double wlen = w.LengthSquare();
// if (wlen <= 0.03) { // ignore walls smaller than 0.15m (15cm)
// return F_wrep;
// }
......@@ -417,18 +433,18 @@ Point GompertzModel::ForceRepWall(Pedestrian* ped, const Line& w, const Point& c
Point pinE; // vorher x1, y1
const JEllipse& E = ped->GetEllipse();
const Point& v = ped->GetV();
if (Distance > J_EPS) {
double min_distance_to_wall = 0.1; // 10 cm
if (Distance > min_distance_to_wall) {
e_iw = dist / Distance;
}
else {
Log->Write("WARNING:\t Gompertz: forceRepWall() ped %d is too near to the wall",ped->GetID());
Log->Write("INFO:\t\t --- take subroom centroid ",ped->GetID());
// Log->Write("WARNING:\t Gompertz: forceRepWall() ped %d is too near to the wall",ped->GetID());
Point new_dist = centroid - ped->GetPos();
new_dist = new_dist/new_dist.Norm();
e_iw = ((inside==true)? new_dist:new_dist*-1);
Distance = EPS;
// Distance = EPS;
// Log->Write("INFO:\t\t --- dist = %f, e= %f %f inside=%d",ped->GetID(), Distance, e_iw.GetX(), e_iw.GetY(), inside);
}
//------------------------- check if others are behind using v0 instead of v
......@@ -438,24 +454,25 @@ Point GompertzModel::ForceRepWall(Pedestrian* ped, const Line& w, const Point& c
// double wallIsBehindv = exp(-exp(-5*tmpv)); //step function: continuous version
double wallIsBehindv = (tmpv<=0)?0:1;
#if DEBUG
if(ped->GetID()==-4){
printf("Distance = %f tmpv=%f\n",Distance, tmpv);
printf("v = %f, %f \n", v.GetX(), v.GetY());
printf("pos = %f, %f \n", ped->GetPos().GetX(), ped->GetPos().GetY());
printf("pt = %f, %f \n", pt.GetX(), pt.GetY());
printf("e_iw = %f, %f\n",e_iw.GetX(), e_iw.GetY());
printf("WallIsBehind = %f (%f)\n",wallIsBehindv,J_EPS);
printf("WallIsBehind = %f (%f)\n",wallIsBehindv,J_EPS);}
#endif
if (wallIsBehindv < J_EPS) { // Wall is behind the direction of motion
if (wallIsBehindv < J_EPS && Distance > min_distance_to_wall) { // Wall is behind the direction of motion
return F_wrep;
}
//------------------------------------------------------------------------
// pt in coordinate system of Ellipse
pinE = pt.TransformToEllipseCoordinates(E.GetCenter(), E.GetCosPhi(), E.GetSinPhi());
pinE = pt.TransformToEllipseCoordinates(E.GetCenter(), E.GetCosPhi(), E.GetSinPhi());
// Punkt auf der Ellipse
r = E.PointOnEllipse(pinE);
Radius = (r - E.GetCenter()).Norm();
// Radius = 0.1;
double radiuss = (r - E.GetCenter()).Norm();
Radius = 0.4;
//-------------------------
const Point& pos = ped->GetPos();
......
......@@ -195,6 +195,7 @@ void Pedestrian::SetV(const Point& v)
_ellipse.SetV(v);
//save the last values for the records
_lastVelocites.push(v);
unsigned int max_size = _recordingTime / _deltaT;
if (_lastVelocites.size() > max_size)
_lastVelocites.pop();
......@@ -448,6 +449,7 @@ double Pedestrian::GetV0Norm() const
}
// orthogonal projection on the stair
//return _ellipse.GetV0()*_building->GetRoom(_roomID)->GetSubRoom(_subRoomID)->GetCosAngleWithHorizontal();
}
// get axis in the walking direction
double Pedestrian::GetLargerAxis() const
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment