modelling wildfire in Surgut

Author: peacekeeper228

include/CellStorage.h

@@ -31,9 +31,6 @@ class CellStorage
 
 public:
 
-    double_t latitudeMin;
-    double_t longtitudeMin;
-
     CellStorage(Math* formula);
     CellStorage(Math* formula, int16_t x_size, int16_t y_size);
     ~CellStorage();

include/Math.h

@@ -16,7 +16,7 @@ class Math
 };
 
 constexpr int ignitionPercentage(){
-    return 150; // in formula there is smth wrong I suppose
+    return 50; // in formula there is smth wrong I suppose
 }
 
 constexpr int throughPercentage(){

include/ProfilingDecorator.h

@@ -39,6 +39,7 @@ class ProfilingDecorator final : public Math
     mutable int positive_counter_;
     mutable int counter_;
     mutable std::unordered_map<metaData, int> results_;
+    mutable std::unordered_map<directions, int> results_directions_;
 
 public:
     ProfilingDecorator(Math *formula);

include/Properties.h

@@ -7,11 +7,11 @@
 
 
 constexpr int getXArea(){
-    return 110;
+    return 1000;
 }
 
 constexpr int getYArea(){
-    return 110;
+    return 1000;
 }
 
 constexpr size_t numberOfSimulations(){
@@ -56,7 +56,8 @@ constexpr int minutesPerIteration(){
 }
 
 constexpr std::string analyzedPolygon(){
-    return "ST_GeomFromText('POLYGON((36.7 55.3, 36.75 55.3, 36.75 55.35, 36.7 55.35, 36.7 55.3))', 4326)";
+    // return "ST_GeomFromText('POLYGON((36.7 55.3, 36.75 55.3, 36.75 55.35, 36.7 55.35, 36.7 55.3))', 4326)";
+    return "ST_GeomFromText('POLYGON((70 60.5, 70.5 60.5, 70.5 61, 70 61, 70 60.5))', 4326)";
 }
 
 constexpr double pi() { return std::atan(1) * 4; }

src/Connection.cpp

@@ -37,12 +37,12 @@ void Connection::setDemToStorage(CellStorage& storage)
         SELECT\
             q.*\
         FROM\
-            (SELECT ST_Rescale(rast, 0.0005, 0.0005) as rast from dem) r,\
+            (SELECT ST_Rescale(rast, 0.00025, 0.00025) as rast from dem_s) r,\
             LATERAL ST_PixelAsCentroids(ST_Clip(r.rast, {}), 1) as q\
         WHERE \
             ST_Intersects(r.rast, {})\
             ) as t\
-    WHERE val > 0;", analyzedPolygon(), analyzedPolygon());
+    WHERE val > 0 and x < {} and y < {};", analyzedPolygon(), analyzedPolygon(),  getXArea(), getYArea());
 
     printf("Query started: %d\n", 1);
     PGresult *res = this->runQuery(formatted_dem_query.c_str());
@@ -67,12 +67,12 @@ void Connection::setBiomassToStorage(CellStorage& storage)
         SELECT\
             q.*\
         FROM\
-            (SELECT ST_Rescale(rast, 0.0005, 0.00025) as rast from biomass) r,\
+            (SELECT ST_Rescale(rast, 0.0005, 0.00025) as rast from biomass_s) r,\
             LATERAL ST_PixelAsCentroids(ST_Clip(r.rast, {}), 1) as q\
         WHERE \
             ST_Intersects(r.rast, {})\
             ) as t\
-    WHERE val > 80 and y < {};", analyzedPolygon(), analyzedPolygon(), getYArea());
+    WHERE val > 80 and x < {} and y < {};", analyzedPolygon(), analyzedPolygon(), getXArea(), getYArea());
     PGresult *res = this->runQuery(formatted_biomass_query.c_str());
     printf("biomass %d\t", PQntuples(res));
     for (int i = 0; i < PQntuples(res); i++)
@@ -87,7 +87,9 @@ void Connection::setBiomassToStorage(CellStorage& storage)
 void Connection::setFireToStorage(CellStorage& storage)
 {
     const std::string formatted_fire_query = std::format("\
-    SELECT brightness, ST_X(geom), ST_Y(geom)\
+    SELECT brightness,\
+        ROUND(ST_Distance(f.geom::geography, ST_SetSRID(ST_MakePoint('70', ST_Y(f.geom)), 4326)::geography)/30),\
+        ROUND(ST_Distance(f.geom::geography, ST_SetSRID(ST_MakePoint(ST_X(f.geom),'60.5'), 4326)::geography)/30)\
     FROM FIRE f\
     WHERE ST_Intersects(f.geom, {})", analyzedPolygon());
     PGresult *res = this->runQuery(formatted_fire_query.c_str());
@@ -96,8 +98,8 @@ void Connection::setFireToStorage(CellStorage& storage)
     {
         double_t k = atof(PQgetvalue(res, i, 1));
         storage.setNewState(cellState::Fire,
-                                  (atof(PQgetvalue(res, i, 1)) - storage.latitudeMin) * 2,
-                                  (atof(PQgetvalue(res, i, 2)) - storage.longtitudeMin) * 2);
+                                  (atoi(PQgetvalue(res, i, 1))),
+                                  (atoi(PQgetvalue(res, i, 2))));
     };
 }
 

src/Math.cpp

@@ -24,8 +24,8 @@ double Math1::calculateWindKoef(const cell *c, directions InvestigatedDirection)
     if (wind == nullptr){
         return double(1);
     }
-    float angleRadians = wind->angleBetweenDirections(InvestigatedDirection) * pi() / 180;
-    double result = std::exp(wind->getWindSpeed() * (cos(cos(angleRadians)) - 1));
+    float angleRadians = wind->angleBetweenDirections(InvestigatedDirection);// * pi() / 180;
+    double result = std::exp(wind->getWindSpeed() * (cos(angleRadians) - 1));
     result_wind_[result]++;
     return result;
 }
@@ -51,7 +51,7 @@ bool Math1::willSpread(const cell *c, directions InvestigatedDirection, int alti
 {
     double fireKoeff = calculateWindKoef(c, InvestigatedDirection);
     double slopeKoeff = calculateGroundSlopeKoef(InvestigatedDirection, altitudeDifference);
-    int result = fireKoeff * slopeKoeff * 100;
+    int result = fireKoeff * slopeKoeff * 100 - 40;
     result_overall[result]++;
     return result + (rand() % 100) > ignitionPercentage();
 }

src/ProfilingDecorator.cpp

@@ -15,6 +15,7 @@ bool ProfilingDecorator::willSpread(const cell *c, directions InvestigatedDirect
     this->counter_++;
     if (result)
     {
+        results_directions_[InvestigatedDirection]++;
         this->positive_counter_++;
     };
     return result;
@@ -46,6 +47,11 @@ ProfilingDecorator::~ProfilingDecorator()
     {
         printf("%d;%d\n", i.first, i.second);
     };
+    printf("Distribution direction\n");
+    for (auto &&i : results_directions_)
+    {
+        printf("%d;%d\n", static_cast<int>(i.first), i.second);
+    };
 
     primary_class_->~Math();
 }

src/Wind.cpp

@@ -36,10 +36,10 @@ double Wind::CalculateWindKoef(directions InvestigatedDirection) const
 int Wind::angleBetweenDirections(directions d) const
 {
     int differenceInDirections = abs(static_cast<int>(d) - static_cast<int>(this->direction_));
-    if (differenceInDirections > 4)
-    {
-        return (8 - differenceInDirections) * numberOfDegreesPerDirection();
-    };
+    // if (differenceInDirections > 4)
+    // {
+    //     return (8 - differenceInDirections) * numberOfDegreesPerDirection();
+    // };
     return differenceInDirections * numberOfDegreesPerDirection();
 }
 

src/main.cpp

@@ -25,11 +25,9 @@ int main()
     Math1 formula = Math1();
     auto formulaDecorator = ProfilingDecorator(&formula);
     CellStorage s = CellStorage(&formulaDecorator);
-    s.latitudeMin = 10.0;
-    s.longtitudeMin = 30.0;
     auto xRange = std::make_pair<int, int>(0, getXArea());
     auto yRange = std::make_pair<int, int>(0, getYArea()/2);
-    auto w = std::make_shared<const Wind>(directions::North, float(1.0));
+    auto w = std::make_shared<const Wind>(directions::SouthWest, float(5.0));
     s.setWindToArea(xRange, yRange, w);
 
     yRange = std::make_pair<int, int>(getYArea()/2, getYArea());