[topography_generator] Added topography generator tool.

16 files changed, 1500 insertions, 0 deletions

diff --git a/topography_generator/CMakeLists.txt b/topography_generator/CMakeLists.txt
new file mode 100644
index 0000000000..6ac22ccfbb
--- /dev/null
+++ b/topography_generator/CMakeLists.txt
@@ -0,0 +1,51 @@
+project(topography_generator_tool)
+
+include_directories(
+  ${OMIM_ROOT}/3party/gflags/src
+)
+
+set(
+  SRC
+  filters_impl.cpp
+  filters_impl.hpp
+  generator.cpp
+  generator.hpp
+  isolines_utils.cpp
+  isolines_utils.hpp
+  main.cpp
+  marching_squares/contours_builder.cpp
+  marching_squares/contours_builder.hpp
+  marching_squares/marching_squares.hpp
+  marching_squares/square.hpp
+  tile_filter.hpp
+  utils/contours.hpp
+  utils/contours_serdes.hpp
+  utils/values_provider.hpp
+)
+
+omim_add_executable(${PROJECT_NAME} ${SRC})
+
+omim_link_libraries(
+  ${PROJECT_NAME}
+  generator
+  routing
+  routing_common
+  traffic
+  storage
+  platform
+  indexer
+  geometry
+  coding
+  base
+  oauthcpp
+  expat
+  jansson
+  protobuf
+  minizip
+  icu
+  stats_client
+  gflags
+  ${LIBZ}
+)
+
+link_qt5_core(${PROJECT_NAME})
diff --git a/topography_generator/filters_impl.cpp b/topography_generator/filters_impl.cpp
new file mode 100644
index 0000000000..6417049707
--- /dev/null
+++ b/topography_generator/filters_impl.cpp
@@ -0,0 +1,24 @@
+#include "topography_generator/filters_impl.hpp"
+
+namespace topography_generator
+{
+std::vector<double> CalculateGaussianLinearKernel(double standardDeviation, double radiusFactor)
+{
+  auto const kernelRadius = static_cast<size_t>(ceil(radiusFactor * standardDeviation));
+  auto const kernelSize = 2 * kernelRadius + 1;
+  std::vector<double> linearKernel(kernelSize, 0);
+
+  double sum = 1.0;
+  linearKernel[kernelRadius] = 1.0;
+  for (int i = 1; i <= kernelRadius; ++i)
+  {
+    double const val = exp(-i * i / (2 * standardDeviation * standardDeviation));
+    linearKernel[kernelRadius - i] = linearKernel[kernelRadius + i] = val;
+    sum += 2.0 * val;
+  }
+  for (auto & val : linearKernel)
+    val /= sum;
+
+  return linearKernel;
+}
+}  // namespace topography_generator
diff --git a/topography_generator/filters_impl.hpp b/topography_generator/filters_impl.hpp
new file mode 100644
index 0000000000..28d210242b
--- /dev/null
+++ b/topography_generator/filters_impl.hpp
@@ -0,0 +1,138 @@
+#pragma once
+
+#include "topography_generator/utils/values_provider.hpp"
+
+namespace topography_generator
+{
+template <typename ValueType>
+void GetExtendedTile(ms::LatLon const & leftBottom, size_t stepsInDegree,
+                     size_t tileSize, size_t tileSizeExtension,
+                     ValuesProvider<ValueType> & valuesProvider,
+                     std::vector<ValueType> & extTileValues)
+{
+  size_t const extendedTileSize = tileSize + 2 * tileSizeExtension;
+  extTileValues.resize(extendedTileSize * extendedTileSize);
+
+  double const step = 1.0 / stepsInDegree;
+  double const offset = step * tileSizeExtension;
+
+  // Store values from North to South.
+  ms::LatLon startPos = ms::LatLon(leftBottom.m_lat + 1.0 + offset,
+                                   leftBottom.m_lon - offset);
+  for (size_t i = 0; i < extendedTileSize; ++i)
+  {
+    for (size_t j = 0; j < extendedTileSize; ++j)
+    {
+      auto const pos = ms::LatLon(startPos.m_lat - i * step,
+                                  startPos.m_lon + j * step);
+      extTileValues[i * extendedTileSize + j] = valuesProvider.GetValue(pos);
+    }
+  }
+}
+
+template <typename ValueType>
+void ProcessWithLinearKernel(std::vector<double> const & kernel, size_t tileSize, size_t tileOffset,
+                             std::vector<ValueType> const & srcValues,
+                             std::vector<ValueType> & dstValues)
+{
+  auto const kernelSize = kernel.size();
+  auto const kernelRadius = kernel.size() / 2;
+  CHECK_LESS_OR_EQUAL(kernelRadius, tileOffset, ());
+  CHECK_GREATER(tileSize, tileOffset * 2, ());
+
+  std::vector<ValueType> tempValues(tileSize, 0);
+
+  for (size_t i = tileOffset; i < tileSize - tileOffset; ++i)
+  {
+    for (size_t j = tileOffset; j < tileSize - tileOffset; ++j)
+    {
+      for (size_t k = 0; k < kernelSize; ++k)
+      {
+        size_t const srcIndex = i * tileSize + j - kernelRadius + k;
+        tempValues[j] += kernel[k] * srcValues[srcIndex];
+      }
+    }
+    for (size_t j = tileOffset; j < tileSize - tileOffset; ++j)
+    {
+      dstValues[i * tileSize + j] = tempValues[j];
+    }
+  }
+
+  for (size_t j = tileOffset; j < tileSize - tileOffset; ++j)
+  {
+    for (size_t i = tileOffset; i < tileSize - tileOffset; ++i)
+    {
+      for (size_t k = 0; k < kernelSize; ++k)
+      {
+        size_t const srcIndex = (i - kernelRadius + k) * tileSize + j;
+        tempValues[i] += kernel[k] * dstValues[srcIndex];
+      }
+    }
+    for (size_t i = tileOffset; i < tileSize - tileOffset; ++i)
+    {
+      dstValues[i * tileSize + j] = tempValues[i];
+    }
+  }
+}
+
+template <typename ValueType>
+void ProcessWithSquareKernel(std::vector<double> const & kernel, size_t kernelSize,
+                             size_t tileSize, size_t tileOffset,
+                             std::vector<ValueType> const & srcValues,
+                             std::vector<ValueType> & dstValues)
+{
+  CHECK_EQUAL(kernelSize * kernelSize, kernel.size(), ());
+  size_t const kernelRadius = kernelSize / 2;
+  CHECK_LESS_OR_EQUAL(kernelRadius, tileOffset, ());
+  CHECK_GREATER(tileSize, tileOffset * 2, ());
+
+  for (size_t i = tileOffset; i < tileSize - tileOffset; ++i)
+  {
+    for (size_t j = tileOffset; j < tileSize - tileOffset; ++j)
+    {
+      size_t const dstIndex = i * tileSize + j;
+      dstValues[dstIndex] = 0;
+      for (size_t ki = 0; ki < kernelSize; ++ki)
+      {
+        for (size_t kj = 0; kj < kernelSize; ++kj)
+        {
+          size_t const srcIndex = (i - kernelRadius + ki) * tileSize + j - kernelRadius + kj;
+          dstValues[dstIndex] += kernel[ki * kernelSize + kj] * srcValues[srcIndex];
+        }
+      }
+    }
+  }
+}
+
+template <typename ValueType>
+void ProcessMedian(size_t kernelRadius, size_t tileSize, size_t tileOffset,
+                   std::vector<ValueType> const & srcValues,
+                   std::vector<ValueType> & dstValues)
+{
+  CHECK_LESS_OR_EQUAL(kernelRadius, tileOffset, ());
+  CHECK_GREATER(tileSize, tileOffset * 2, ());
+
+  size_t const kernelSize = kernelRadius * 2 + 1;
+  std::vector<ValueType> kernel(kernelSize * kernelSize);
+  for (size_t i = tileOffset; i < tileSize - tileOffset; ++i)
+  {
+    for (size_t j = tileOffset; j < tileSize - tileOffset; ++j)
+    {
+      size_t const startI = i - kernelRadius;
+      size_t const startJ = j - kernelRadius;
+      for (size_t ki = 0; ki < kernelSize; ++ki)
+      {
+        for (size_t kj = 0; kj < kernelSize; ++kj)
+        {
+          size_t const srcIndex = (startI + ki) * tileSize + startJ + kj;
+          kernel[ki * kernelSize + kj] = srcValues[srcIndex];
+        }
+      }
+      std::sort(kernel.begin(), kernel.end());
+      dstValues[i * tileSize + j] = kernel[kernelRadius];
+    }
+  }
+}
+
+std::vector<double> CalculateGaussianLinearKernel(double standardDeviation, double radiusFactor);
+}  // namespace topography_generator
diff --git a/topography_generator/generator.cpp b/topography_generator/generator.cpp
new file mode 100644
index 0000000000..6a997aa627
--- /dev/null
+++ b/topography_generator/generator.cpp
@@ -0,0 +1,282 @@
+#include "topography_generator/generator.hpp"
+#include "topography_generator/isolines_utils.hpp"
+#include "topography_generator/marching_squares/marching_squares.hpp"
+#include "topography_generator/utils/contours_serdes.hpp"
+
+#include "platform/platform.hpp"
+
+#include "generator/srtm_parser.hpp"
+
+#include "geometry/mercator.hpp"
+
+#include <vector>
+
+namespace topography_generator
+{
+size_t constexpr kArcSecondsInDegree = 60 * 60;
+
+class SrtmProvider : public ValuesProvider<Altitude>
+{
+public:
+  explicit SrtmProvider(std::string const & srtmDir):
+    m_srtmManager(srtmDir)
+  {}
+
+  Altitude GetValue(ms::LatLon const & pos) override
+  {
+    return m_srtmManager.GetHeight(pos);
+  }
+
+  Altitude GetInvalidValue() const override
+  {
+    return kInvalidAltitude;
+  }
+
+private:
+  generator::SrtmTileManager m_srtmManager;
+};
+
+class RawSrtmTile : public ValuesProvider<Altitude>
+{
+public:
+  RawSrtmTile(std::vector<Altitude> const & values,
+              int leftLon, int bottomLat)
+    : m_values(values)
+    , m_leftLon(leftLon)
+    , m_bottomLat(bottomLat)
+  {}
+
+  Altitude GetValue(ms::LatLon const & pos) override
+  {
+    // TODO: assert
+    CHECK_EQUAL(floor(pos.m_lat), m_bottomLat, ());
+    CHECK_EQUAL(floor(pos.m_lon), m_leftLon, ());
+
+    double ln = pos.m_lon - m_leftLon;
+    double lt = pos.m_lat - m_bottomLat;
+    lt = 1 - lt;  // from North to South
+
+    size_t const row = kArcSecondsInDegree * lt + 0.5;
+    size_t const col = kArcSecondsInDegree * ln + 0.5;
+
+    size_t const ix = row * (kArcSecondsInDegree + 1) + col;
+    return ix < m_values.size() ? m_values[ix] : kInvalidAltitude;
+  }
+
+  Altitude GetInvalidValue() const override
+  {
+    return kInvalidAltitude;
+  }
+
+private:
+  std::vector<Altitude> const & m_values;
+  int m_leftLon;
+  int m_bottomLat;
+};
+
+class TileIsolinesTask : public threads::IRoutine
+{
+public:
+  TileIsolinesTask(int left, int bottom, int right, int top,
+                   std::string const & srtmDir, TileIsolinesParams const & params)
+    : m_left(left)
+    , m_bottom(bottom)
+    , m_right(right)
+    , m_top(top)
+    , m_srtmProvider(srtmDir)
+    , m_params(params)
+  {}
+
+  void Do() override
+  {
+    for (int lat = m_bottom; lat <= m_top; ++lat)
+    {
+      for (int lon = m_left; lon <= m_right; ++lon)
+        ProcessTile(lat, lon);
+    }
+  }
+
+private:
+  void ProcessTile(int lat, int lon)
+  {
+    auto const tileName = generator::SrtmTile::GetBase(ms::LatLon(lat, lon));
+    LOG(LINFO, ("Begin generating isolines for tile", tileName));
+    ms::LatLon const leftBottom = ms::LatLon(lat, lon);
+    ms::LatLon const rightTop = ms::LatLon(lat + 1.0, lon + 1.0);
+    double const squaresStep = 1.0 / (kArcSecondsInDegree) * m_params.m_latLonStepFactor;
+    Contours<Altitude> contours;
+    if (lat >= 60)
+    {
+      // Filter tiles converted from ASTER, cause they are noisy enough.
+      std::vector<Altitude> filteredValues = FilterTile(
+        m_params.m_filters, leftBottom, kArcSecondsInDegree,
+        kArcSecondsInDegree + 1, m_srtmProvider);
+      RawSrtmTile filteredProvider(filteredValues, lon, lat);
+
+      MarchingSquares<Altitude> squares(leftBottom, rightTop,
+                                        squaresStep, m_params.m_alitudesStep,
+                                        filteredProvider);
+      squares.GenerateContours(contours);
+    }
+    else
+    {
+      MarchingSquares<Altitude> squares(leftBottom, rightTop,
+                                        squaresStep, m_params.m_alitudesStep,
+                                        m_srtmProvider);
+      squares.GenerateContours(contours);
+    }
+
+    SaveContrours(GetIsolinesFilePath(lat, lon, m_params.m_outputDir), std::move(contours));
+    LOG(LINFO, ("End generating isolines for tile", tileName));
+  }
+
+  int m_left;
+  int m_bottom;
+  int m_right;
+  int m_top;
+  SrtmProvider m_srtmProvider;
+  TileIsolinesParams const & m_params;
+};
+
+Generator::Generator(std::string const & srtmPath, size_t threadsCount, size_t maxCachedTilesPerThread)
+  : m_threadsCount(threadsCount)
+  , m_maxCachedTilesPerThread(maxCachedTilesPerThread)
+  , m_srtmPath(srtmPath)
+{
+  m_infoGetter = storage::CountryInfoReader::CreateCountryInfoReader(GetPlatform());
+  CHECK(m_infoGetter, ());
+  m_infoReader = dynamic_cast<storage::CountryInfoReader *>(m_infoGetter.get());
+
+  m_threadsPool = std::make_unique<base::thread_pool::routine::ThreadPool>(
+    threadsCount, std::bind(&Generator::OnTaskFinished, this, std::placeholders::_1));
+}
+
+Generator::~Generator()
+{
+  m_threadsPool->Stop();
+}
+
+void Generator::GenerateIsolines(int left, int bottom, int right, int top,
+                                 TileIsolinesParams const & params)
+{
+  std::vector<std::unique_ptr<TileIsolinesTask>> tasks;
+
+  CHECK_GREATER_OR_EQUAL(right, left, ());
+  CHECK_GREATER_OR_EQUAL(top, bottom, ());
+
+  int tilesRowPerTask = top - bottom + 1;
+  int tilesColPerTask = right - left + 1;
+
+  if (tilesRowPerTask * tilesColPerTask <= m_threadsCount)
+  {
+    tilesRowPerTask = 1;
+    tilesColPerTask = 1;
+  }
+  else
+  {
+    while (tilesRowPerTask * tilesColPerTask > m_maxCachedTilesPerThread)
+    {
+      if (tilesRowPerTask > tilesColPerTask)
+        tilesRowPerTask = (tilesRowPerTask + 1) / 2;
+      else
+        tilesColPerTask = (tilesColPerTask + 1) / 2;
+    }
+  }
+
+  for (int lat = bottom; lat <= top; lat += tilesRowPerTask)
+  {
+    int const topLat = std::min(lat + tilesRowPerTask - 1, top);
+    for (int lon = left; lon <= right; lon += tilesColPerTask)
+    {
+      int const rightLon = std::min(lon + tilesColPerTask - 1, right);
+      tasks.emplace_back(new TileIsolinesTask(lon, lat, rightLon, topLat, m_srtmPath, params));
+    }
+  }
+
+  {
+    std::lock_guard<std::mutex> lock(m_tasksMutex);
+    CHECK(m_activeTasksCount == 0, ());
+    m_activeTasksCount = tasks.size();
+  }
+
+  for (auto & task : tasks)
+    m_threadsPool->PushBack(task.get());
+
+  std::unique_lock<std::mutex> lock(m_tasksMutex);
+  m_tasksReadyCondition.wait(lock, [this] { return m_activeTasksCount == 0; });
+}
+
+void Generator::OnTaskFinished(threads::IRoutine * task)
+{
+  {
+    std::lock_guard<std::mutex> lock(m_tasksMutex);
+    CHECK(m_activeTasksCount > 0, ());
+    --m_activeTasksCount;
+    if (m_activeTasksCount == 0)
+      m_tasksReadyCondition.notify_one();
+  }
+}
+
+void Generator::GetCountryRegions(std::string const & countryId, m2::RectD & countryRect,
+                                  std::vector<m2::RegionD> & countryRegions)
+{
+  countryRect = m_infoReader->GetLimitRectForLeaf(countryId);
+
+  size_t id;
+  for (id = 0; id < m_infoReader->GetCountries().size(); ++id)
+  {
+    if (m_infoReader->GetCountries().at(id).m_countryId == countryId)
+      break;
+  }
+  CHECK_LESS(id, m_infoReader->GetCountries().size(), ());
+
+  m_infoReader->LoadRegionsFromDisk(id, countryRegions);
+}
+
+void Generator::PackIsolinesForCountry(std::string const & countryId,
+                                       std::string const & isolinesPath,
+                                       std::string const & outDir,
+                                       CountryIsolinesParams const & params)
+{
+  m2::RectD countryRect;
+  std::vector<m2::RegionD> countryRegions;
+  GetCountryRegions(countryId, countryRect, countryRegions);
+
+  auto const leftBottom = mercator::ToLatLon(countryRect.LeftBottom());
+  auto const rightTop = mercator::ToLatLon(countryRect.RightTop());
+
+  auto const left = static_cast<int>(floor(leftBottom.m_lon));
+  auto const bottom = static_cast<int>(floor(leftBottom.m_lat));
+  auto const right = static_cast<int>(floor(rightTop.m_lon));
+  auto const top = static_cast<int>(floor(rightTop.m_lat));
+
+  Contours<Altitude> countryIsolines;
+  countryIsolines.m_minValue = std::numeric_limits<Altitude>::max();
+  countryIsolines.m_maxValue = std::numeric_limits<Altitude>::min();
+
+  for (int lat = bottom; lat <= top; ++lat)
+  {
+    for (int lon = left; lon <= right; ++lon)
+    {
+      Contours<Altitude> isolines;
+      if (!LoadContours(GetIsolinesFilePath(lat, lon, isolinesPath), isolines))
+        continue;
+
+      if (params.m_simplificationZoom > 0)
+        SimplifyContours(params.m_simplificationZoom, isolines);
+      CropContours(countryRect, countryRegions, params.m_maxIsolineLenght, isolines);
+
+      countryIsolines.m_minValue = std::min(isolines.m_minValue, countryIsolines.m_minValue);
+      countryIsolines.m_maxValue = std::max(isolines.m_maxValue, countryIsolines.m_maxValue);
+      countryIsolines.m_valueStep = isolines.m_valueStep;
+      for (auto & levelIsolines : isolines.m_contours)
+      {
+        auto & dst = countryIsolines.m_contours[levelIsolines.first];
+        dst.insert(dst.end(), levelIsolines.second.begin(), levelIsolines.second.end());
+      }
+    }
+  }
+
+  SaveContrours(GetIsolinesFilePath(countryId, outDir), std::move(countryIsolines));
+}
+}  // namespace topography_generator
diff --git a/topography_generator/generator.hpp b/topography_generator/generator.hpp
new file mode 100644
index 0000000000..daaa219505
--- /dev/null
+++ b/topography_generator/generator.hpp
@@ -0,0 +1,61 @@
+#pragma once
+
+#include "topography_generator/isolines_utils.hpp"
+#include "topography_generator/tile_filter.hpp"
+
+#include "storage/country_info_getter.hpp"
+
+#include "geometry/rect2d.hpp"
+#include "geometry/region2d.hpp"
+
+#include "base/thread_pool.hpp"
+
+#include <condition_variable>
+#include <string>
+#include <memory>
+#include <mutex>
+
+namespace topography_generator
+{
+struct TileIsolinesParams
+{
+  Altitude m_alitudesStep = 10;
+  size_t m_latLonStepFactor = 1;
+  FiltersSequence<Altitude> m_filters;
+  std::string m_outputDir;
+};
+
+struct CountryIsolinesParams
+{
+  size_t m_maxIsolineLenght = 1000;
+  int m_simplificationZoom = 17;
+};
+
+class Generator
+{
+public:
+  Generator(std::string const & srtmPath, size_t threadsCount, size_t maxCachedTilesPerThread);
+  ~Generator();
+
+  void GenerateIsolines(int left, int bottom, int right, int top,
+                        TileIsolinesParams const & params);
+  void PackIsolinesForCountry(std::string const & countryId, std::string const & isolinesPath,
+                              std::string const & outDir, CountryIsolinesParams const & params);
+
+private:
+  void OnTaskFinished(threads::IRoutine * task);
+  void GetCountryRegions(std::string const & countryId, m2::RectD & countryRect,
+                         std::vector<m2::RegionD> & countryRegions);
+
+  std::unique_ptr<storage::CountryInfoGetter> m_infoGetter;
+  storage::CountryInfoReader * m_infoReader = nullptr;
+
+  std::unique_ptr<base::thread_pool::routine::ThreadPool> m_threadsPool;
+  size_t m_threadsCount;
+  size_t m_maxCachedTilesPerThread;
+  std::string m_srtmPath;
+  std::mutex m_tasksMutex;
+  std::condition_variable m_tasksReadyCondition;
+  size_t m_activeTasksCount = 0;
+};
+}  // namespace topography_generator
diff --git a/topography_generator/isolines_utils.cpp b/topography_generator/isolines_utils.cpp
new file mode 100644
index 0000000000..e07fa4a70d
--- /dev/null
+++ b/topography_generator/isolines_utils.cpp
@@ -0,0 +1,21 @@
+#include "topography_generator/isolines_utils.hpp"
+
+#include "generator/srtm_parser.hpp"
+
+#include "base/file_name_utils.hpp"
+
+namespace topography_generator
+{
+std::string const kIsolinesExt = ".isolines";
+
+std::string GetIsolinesFilePath(int lat, int lon, std::string const & dir)
+{
+  auto const fileName = generator::SrtmTile::GetBase(ms::LatLon(lat, lon));
+  return base::JoinPath(dir, fileName + kIsolinesExt);
+}
+
+std::string GetIsolinesFilePath(std::string const & countryId, std::string const & dir)
+{
+  return base::JoinPath(dir, countryId + kIsolinesExt);
+}
+}  // namespace topography_generator
diff --git a/topography_generator/isolines_utils.hpp b/topography_generator/isolines_utils.hpp
new file mode 100644
index 0000000000..d2c9720be8
--- /dev/null
+++ b/topography_generator/isolines_utils.hpp
@@ -0,0 +1,14 @@
+#pragma once
+
+#include "indexer/feature_altitude.hpp"
+
+#include <string>
+
+namespace topography_generator
+{
+using Altitude = feature::TAltitude;
+Altitude constexpr kInvalidAltitude = feature::kInvalidAltitude;
+
+std::string GetIsolinesFilePath(int lat, int lon, std::string const & dir);
+std::string GetIsolinesFilePath(std::string const & countryId, std::string const & dir);
+}  // namespace topography_generator
diff --git a/topography_generator/main.cpp b/topography_generator/main.cpp
new file mode 100644
index 0000000000..bcf983732c
--- /dev/null
+++ b/topography_generator/main.cpp
@@ -0,0 +1,87 @@
+#include "generator.hpp"
+
+#include "base/assert.hpp"
+#include "base/logging.hpp"
+
+#include "3party/gflags/src/gflags/gflags.h"
+
+
+DEFINE_string(out_dir, "/Users/daravolvenkova/isolines/", "Path to output directory.");
+DEFINE_string(countryId, "",//"France_Rhone-Alpes_Haute-Savoie",
+              "Isolines packing mode. Pack isolines for countryId.");
+DEFINE_string(isolines_path, "/Users/daravolvenkova/isolines/",
+              "Isolines packing mode. Path to the directory with isolines tiles.");
+DEFINE_uint64(simpl_zoom, 17, "Isolines packing mode. Isolines simplification zoom.");
+DEFINE_uint64(max_length, 1000, "Isolines packing mode. Isolines max length.");
+DEFINE_string(srtm_path, "/Users/daravolvenkova/srtm/2000.02.11/",
+              "Isolines generating mode. Path to srtm directory.");
+DEFINE_int32(left, 5, "Isolines generating mode. Left longitude of tiles rect.");
+DEFINE_int32(right, 6, "Isolines generating mode. Right longitude of tiles rect.");
+DEFINE_int32(bottom, 45, "Isolines generating mode. Bottom latitude of tiles rect.");
+DEFINE_int32(top, 46, "Isolines generating mode. Top latitude of tiles rect.");
+DEFINE_uint64(isolines_step, 10, "Isolines generating mode. Isolines step in meters.");
+DEFINE_uint64(latlon_step_factor, 1, "Isolines generating mode. Lat/lon step factor.");
+DEFINE_double(gaussian_st_dev, 2.0, "Isolines generating mode. Gaussian filter standard deviation.");
+DEFINE_double(gaussian_r_factor, 1.0, "Isolines generating mode. Gaussian filter radius factor.");
+DEFINE_uint64(median_r, 1, "Isolines generating mode. Median filter radius.");
+DEFINE_uint64(threads, 4, "Number of threads.");
+DEFINE_uint64(tiles_per_thread, 9, "Max cached tiles per thread");
+
+
+int main(int argc, char ** argv)
+{
+  google::ParseCommandLineFlags(&argc, &argv, true);
+
+  if (FLAGS_out_dir.empty())
+  {
+    LOG(LINFO, ("out_dir must be set."));
+    return 1;
+  }
+
+  topography_generator::Generator generator(FLAGS_srtm_path, FLAGS_threads,
+                                            FLAGS_tiles_per_thread);
+  if (!FLAGS_countryId.empty())
+  {
+    if (FLAGS_isolines_path.empty())
+    {
+      LOG(LINFO, ("isolines_path must be set."));
+      return 1;
+    }
+
+    topography_generator::CountryIsolinesParams params;
+    params.m_simplificationZoom = static_cast<int>(FLAGS_simpl_zoom);
+    params.m_maxIsolineLenght = FLAGS_max_length;
+
+    generator.PackIsolinesForCountry(FLAGS_countryId, FLAGS_isolines_path, FLAGS_out_dir, params);
+    return 0;
+  }
+
+  if (FLAGS_srtm_path.empty())
+  {
+    LOG(LINFO, ("srtm_path must be set."));
+    return 1;
+  }
+
+  if (FLAGS_right < FLAGS_left || FLAGS_top < FLAGS_bottom)
+  {
+    LOG(LINFO, ("Invalid tiles rect."));
+    return 1;
+  }
+
+  topography_generator::TileIsolinesParams params;
+  if (FLAGS_median_r > 0)
+    params.m_filters.emplace_back(new topography_generator::MedianFilter<topography_generator::Altitude>(FLAGS_median_r));
+  if (FLAGS_gaussian_st_dev > 0.0)
+  {
+    params.m_filters.emplace_back(new topography_generator::GaussianFilter<topography_generator::Altitude>(
+      FLAGS_gaussian_st_dev, FLAGS_gaussian_r_factor));
+  }
+
+  params.m_outputDir = FLAGS_out_dir;
+  params.m_alitudesStep = FLAGS_isolines_step;
+  params.m_latLonStepFactor = FLAGS_latlon_step_factor;
+
+  generator.GenerateIsolines(FLAGS_left, FLAGS_bottom, FLAGS_right, FLAGS_top, params);
+
+  return 0;
+}
diff --git a/topography_generator/marching_squares/contours_builder.cpp b/topography_generator/marching_squares/contours_builder.cpp
new file mode 100644
index 0000000000..4c5156894e
--- /dev/null
+++ b/topography_generator/marching_squares/contours_builder.cpp
@@ -0,0 +1,102 @@
+#include "topography_generator/marching_squares/contours_builder.hpp"
+
+#include "geometry/mercator.hpp"
+
+namespace topography_generator
+{
+double constexpr kEps = 1e-7;
+
+ContoursBuilder::ContoursBuilder(size_t levelsCount)
+  : m_levelsCount(levelsCount)
+{
+  m_finalizedContours.resize(m_levelsCount);
+  m_activeContours.resize(m_levelsCount);
+}
+
+void ContoursBuilder::AddSegment(size_t levelInd, ms::LatLon const & beginPos, ms::LatLon const & endPos)
+{
+  if (beginPos.EqualDxDy(endPos, kEps))
+    return;
+
+  CHECK_LESS(levelInd, m_levelsCount, ());
+
+  auto contourItBefore = FindContourWithEndPoint(levelInd, beginPos);
+  auto contourItAfter = FindContourWithStartPoint(levelInd, endPos);
+  auto const connectStart = contourItBefore != m_activeContours[levelInd].end();
+  auto const connectEnd = contourItAfter != m_activeContours[levelInd].end();
+
+  if (connectStart && connectEnd && contourItBefore != contourItAfter)
+  {
+    contourItBefore->m_countour.insert(contourItBefore->m_countour.end(),
+      contourItAfter->m_countour.begin(), contourItAfter->m_countour.end());
+    contourItBefore->m_active = true;
+    m_activeContours[levelInd].erase(contourItAfter);
+  }
+  else if (connectStart)
+  {
+    contourItBefore->m_countour.push_back(endPos);
+    contourItBefore->m_active = true;
+  }
+  else if (connectEnd)
+  {
+    contourItAfter->m_countour.push_front(beginPos);
+    contourItBefore->m_active = true;
+  }
+  else
+  {
+    m_activeContours[levelInd].emplace_back(ContourRaw({beginPos, endPos}));
+  }
+}
+
+void ContoursBuilder::BeginLine()
+{
+  for (auto & contoursList : m_activeContours)
+  {
+    for (auto & activeContour : contoursList)
+      activeContour.m_active = false;
+  }
+}
+
+void ContoursBuilder::EndLine(bool finalLine)
+{
+  for (size_t levelInd = 0; levelInd < m_activeContours.size(); ++levelInd)
+  {
+    auto & contours = m_activeContours[levelInd];
+    auto it = contours.begin();
+    while (it != contours.end())
+    {
+      if (!it->m_active || finalLine)
+      {
+        m_finalizedContours[levelInd].push_back(std::move(it->m_countour));
+        it = contours.erase(it);
+      }
+      else
+      {
+        ++it;
+      }
+    }
+  }
+}
+
+ContoursBuilder::ActiveContourIter ContoursBuilder::FindContourWithStartPoint(size_t levelInd, ms::LatLon const & pos)
+{
+  auto & contours = m_activeContours[levelInd];
+  for (auto it = contours.begin(); it != contours.end(); ++it)
+  {
+    if (it->m_countour.front().EqualDxDy(pos, kEps))
+      return it;
+  }
+  return contours.end();
+}
+
+ContoursBuilder::ActiveContourIter ContoursBuilder::FindContourWithEndPoint(size_t levelInd, ms::LatLon const & pos)
+{
+  auto & contours = m_activeContours[levelInd];
+  for (auto it = contours.begin(); it != contours.end(); ++it)
+  {
+    if (it->m_countour.back().EqualDxDy(pos, kEps))
+      return it;
+  }
+  return contours.end();
+}
+}  // namespace topography_generator
diff --git a/topography_generator/marching_squares/contours_builder.hpp b/topography_generator/marching_squares/contours_builder.hpp
new file mode 100644
index 0000000000..4f2a3ea8c0
--- /dev/null
+++ b/topography_generator/marching_squares/contours_builder.hpp
@@ -0,0 +1,67 @@
+#pragma once
+
+#include "topography_generator/utils/contours.hpp"
+
+#include "geometry/latlon.hpp"
+
+#include <deque>
+#include <list>
+#include <vector>
+#include <unordered_map>
+
+namespace topography_generator
+{
+class ContoursBuilder
+{
+public:
+  explicit ContoursBuilder(size_t levelsCount);
+
+  void AddSegment(size_t levelInd, ms::LatLon const & beginPos, ms::LatLon const & endPos);
+  void BeginLine();
+  void EndLine(bool finalLine);
+
+  template <typename ValueType>
+  void GetContours(ValueType minValue, ValueType valueStep,
+                   std::unordered_map<ValueType, std::vector<Contour>> & contours)
+  {
+    contours.clear();
+    for (size_t i = 0; i < m_finalizedContours.size(); ++i)
+    {
+      auto const levelValue = minValue + i * valueStep;
+      auto const & contoursList = m_finalizedContours[i];
+      for (auto const & contour : contoursList)
+      {
+        Contour contourMerc;
+        contourMerc.reserve(contour.size());
+        for (auto const & pt : contour)
+          contourMerc.push_back(mercator::FromLatLon(pt));
+        contours[levelValue].emplace_back(std::move(contourMerc));
+      }
+    }
+  }
+
+private:
+  using ContourRaw = std::deque<ms::LatLon>;
+  using ContoursList = std::list<ContourRaw>;
+
+  struct ActiveContour
+  {
+    explicit ActiveContour(ContourRaw && isoline)
+      : m_countour(std::move(isoline))
+    {}
+
+    ContourRaw m_countour;
+    bool m_active = true;
+  };
+  using ActiveContoursList = std::list<ActiveContour>;
+  using ActiveContourIter = ActiveContoursList::iterator;
+
+  ActiveContourIter FindContourWithStartPoint(size_t levelInd, ms::LatLon const & pos);
+  ActiveContourIter FindContourWithEndPoint(size_t levelInd, ms::LatLon const & pos);
+
+  size_t const m_levelsCount;
+
+  std::vector<ContoursList> m_finalizedContours;
+  std::vector<ActiveContoursList> m_activeContours;
+};
+}  // namespace topography_generator
diff --git a/topography_generator/marching_squares/marching_squares.hpp b/topography_generator/marching_squares/marching_squares.hpp
new file mode 100644
index 0000000000..a57bd55856
--- /dev/null
+++ b/topography_generator/marching_squares/marching_squares.hpp
@@ -0,0 +1,110 @@
+#pragma once
+
+#include "topography_generator/marching_squares/contours_builder.hpp"
+#include "topography_generator/marching_squares/square.hpp"
+#include "topography_generator/utils/contours.hpp"
+
+#include "base/logging.hpp"
+
+namespace topography_generator
+{
+template <typename ValueType>
+class MarchingSquares
+{
+public:
+  MarchingSquares(ms::LatLon const & leftBottom, ms::LatLon const & rightTop,
+                  double step, ValueType valueStep, ValuesProvider<ValueType> & valuesProvider)
+    : m_leftBottom(leftBottom)
+    , m_rightTop(rightTop)
+    , m_step(step)
+    , m_valueStep(valueStep)
+    , m_valuesProvider(valuesProvider)
+  {
+    CHECK_GREATER(m_rightTop.m_lon, m_leftBottom.m_lon, ());
+    CHECK_GREATER(m_rightTop.m_lat, m_leftBottom.m_lat, ());
+
+    m_stepsCountLon = static_cast<size_t>((m_rightTop.m_lon - m_leftBottom.m_lon) / step);
+    m_stepsCountLat = static_cast<size_t>((m_rightTop.m_lat - m_leftBottom.m_lat) / step);
+
+    CHECK_GREATER(m_stepsCountLon, 0, ());
+    CHECK_GREATER(m_stepsCountLat, 0, ());
+  }
+
+  void GenerateContours(Contours<ValueType> & result)
+  {
+    ScanValuesInRect(result.m_minValue, result.m_maxValue, result.m_invalidValuesCount);
+    result.m_valueStep = m_valueStep;
+
+    auto const levelsCount = static_cast<size_t>(result.m_maxValue - result.m_minValue) / m_valueStep;
+    if (levelsCount == 0)
+    {
+      LOG(LINFO, ("Contours can't be generated: min and max values are equal:", result.m_minValue));
+      return;
+    }
+
+    ContoursBuilder contoursBuilder(levelsCount);
+
+    for (size_t i = 0; i < m_stepsCountLat - 1; ++i)
+    {
+      contoursBuilder.BeginLine();
+      for (size_t j = 0; j < m_stepsCountLon - 1; ++j)
+      {
+        auto const pos = ms::LatLon(m_leftBottom.m_lat + m_step * i, m_leftBottom.m_lon + m_step * j);
+        Square<ValueType> square(pos, m_step, result.m_minValue, m_valueStep, m_valuesProvider);
+        square.GenerateSegments(contoursBuilder);
+      }
+      auto const isLastLine = i == m_stepsCountLat - 2;
+      contoursBuilder.EndLine(isLastLine);
+    }
+
+    contoursBuilder.GetContours(result.m_minValue, result.m_valueStep, result.m_contours);
+  }
+
+private:
+  void ScanValuesInRect(ValueType & minValue, ValueType & maxValue, size_t & invalidValuesCount) const
+  {
+    minValue = maxValue = m_valuesProvider.GetValue(m_leftBottom);
+    invalidValuesCount = 0;
+
+    for (size_t i = 0; i < m_stepsCountLat; ++i)
+    {
+      for (size_t j = 0; j < m_stepsCountLon; ++j)
+      {
+        auto const pos = ms::LatLon(m_leftBottom.m_lat + m_step * i,
+                                    m_leftBottom.m_lon + m_step * j);
+        auto const value = m_valuesProvider.GetValue(pos);
+        if (value == m_valuesProvider.GetInvalidValue())
+        {
+          ++invalidValuesCount;
+          continue;
+        }
+        if (value < minValue)
+          minValue = value;
+        if (value > maxValue)
+          maxValue = value;
+      }
+    }
+
+    if (minValue > 0)
+      minValue = m_valueStep * ((minValue + m_valueStep - 1) / m_valueStep);
+    else
+      minValue = m_valueStep * (minValue / m_valueStep);
+
+    if (maxValue > 0)
+      maxValue = m_valueStep * ((maxValue + m_valueStep) / m_valueStep);
+    else
+      maxValue = m_valueStep * (maxValue / m_valueStep);
+
+    CHECK_GREATER_OR_EQUAL(maxValue, minValue, ());
+  }
+
+  ms::LatLon const m_leftBottom;
+  ms::LatLon const m_rightTop;
+  double const m_step;
+  ValueType const m_valueStep;
+  ValuesProvider<ValueType> & m_valuesProvider;
+
+  size_t m_stepsCountLon;
+  size_t m_stepsCountLat;
+};
+}  // namespace topography_generator
diff --git a/topography_generator/marching_squares/square.hpp b/topography_generator/marching_squares/square.hpp
new file mode 100644
index 0000000000..9500a66ab1
--- /dev/null
+++ b/topography_generator/marching_squares/square.hpp
@@ -0,0 +1,207 @@
+#pragma once
+
+#include "topography_generator/marching_squares/contours_builder.hpp"
+#include "topography_generator/utils/values_provider.hpp"
+
+namespace topography_generator
+{
+template <typename ValueType>
+class Square
+{
+public:
+  Square(ms::LatLon const & leftBottom, double size,
+         ValueType minValue, ValueType valueStep,
+         ValuesProvider<ValueType> & valuesProvider)
+    : m_minValue(minValue)
+    , m_valueStep(valueStep)
+    , m_left(leftBottom.m_lon)
+    , m_right(leftBottom.m_lon + size)
+    , m_bottom(leftBottom.m_lat)
+    , m_top(leftBottom.m_lat + size)
+  {
+    static_assert(std::is_integral<ValueType>::value, "Only integral types are supported.");
+
+    m_valueLB = GetValue(leftBottom, valuesProvider);
+    m_valueLT = GetValue(ms::LatLon(m_top, m_left), valuesProvider);
+    m_valueRT = GetValue(ms::LatLon(m_top, m_right), valuesProvider);
+    m_valueRB = GetValue(ms::LatLon(m_bottom, m_right), valuesProvider);
+  }
+
+  void GenerateSegments(ContoursBuilder & builder)
+  {
+    ValueType minAlt = std::min(m_valueLB, std::min(m_valueLT, std::min(m_valueRT, m_valueRB)));
+    ValueType maxAlt = std::max(m_valueLB, std::max(m_valueLT, std::max(m_valueRT, m_valueRB)));
+
+    if (minAlt > 0)
+      minAlt = m_valueStep * ((minAlt + m_valueStep - 1) / m_valueStep);
+    else
+      minAlt = m_valueStep * (minAlt / m_valueStep);
+    if (maxAlt > 0)
+      maxAlt = m_valueStep * ((maxAlt + m_valueStep) / m_valueStep);
+    else
+      maxAlt = m_valueStep * (maxAlt / m_valueStep);
+
+    CHECK_GREATER_OR_EQUAL(minAlt, m_minValue, ());
+
+    for (auto alt = minAlt; alt < maxAlt; alt += m_valueStep)
+      AddSegments(alt, (alt - m_minValue) / m_valueStep, builder);
+  }
+
+//private:
+  enum class Rib
+  {
+    None,
+    Left,
+    Top,
+    Right,
+    Bottom,
+    Unclear,
+  };
+
+  ValueType GetValue(ms::LatLon const & pos, ValuesProvider<ValueType> & valuesProvider) const
+  {
+    ValueType val = valuesProvider.GetValue(pos);
+    if (val != valuesProvider.GetInvalidValue() && (val % m_valueStep == 0))
+      return val + 1;
+    return val;
+  }
+
+  void AddSegments(ValueType val, uint16_t ind, ContoursBuilder & builder)
+  {
+    // Segment is a vector directed so that higher values is on the right.
+    std::pair<Rib, Rib> intersectedRibs[] =
+      {
+        {Rib::None, Rib::None},       // 0000
+        {Rib::Left, Rib::Bottom},     // 0001
+        {Rib::Top, Rib::Left},        // 0010
+        {Rib::Top, Rib::Bottom},      // 0011
+        {Rib::Right, Rib::Top},       // 0100
+        {Rib::Unclear, Rib::Unclear}, // 0101
+        {Rib::Right, Rib::Left},      // 0110
+        {Rib::Right, Rib::Bottom},    // 0111
+        {Rib::Bottom, Rib::Right},    // 1000
+        {Rib::Left, Rib::Right},      // 1001
+        {Rib::Unclear, Rib::Unclear}, // 1010
+        {Rib::Top, Rib::Right},       // 1011
+        {Rib::Bottom, Rib::Top},      // 1100
+        {Rib::Left, Rib::Top},        // 1101
+        {Rib::Bottom, Rib::Left},     // 1110
+        {Rib::None, Rib::None},       // 1111
+      };
+
+    uint8_t const pattern =  (m_valueLB > val ? 1u : 0u) | ((m_valueLT > val ? 1u : 0u) << 1u) |
+      ((m_valueRT > val ? 1u : 0u) << 2u) | ((m_valueRB > val ? 1u : 0u) << 3u);
+
+    auto ribs = intersectedRibs[pattern];
+
+    if (ribs.first == Rib::None)
+      return;
+
+    if (ribs.first != Rib::Unclear)
+    {
+      builder.AddSegment(ind, InterpolatePoint(ribs.first, val), InterpolatePoint(ribs.second, val));
+    }
+    else
+    {
+      auto const leftPos = InterpolatePoint(Rib::Left, val);
+      auto const rightPos = InterpolatePoint(Rib::Right, val);
+      auto const bottomPos = InterpolatePoint(Rib::Bottom, val);
+      auto const topPos = InterpolatePoint(Rib::Top, val);
+
+      ValueType const avVal = (m_valueLB + m_valueLT + m_valueRT + m_valueRB) / 4;
+      if (avVal > val)
+      {
+        if (m_valueLB > val)
+        {
+          builder.AddSegment(ind, leftPos, topPos);
+          builder.AddSegment(ind, rightPos, bottomPos);
+        }
+        else
+        {
+          builder.AddSegment(ind, bottomPos, leftPos);
+          builder.AddSegment(ind, topPos, rightPos);
+        }
+      }
+      else
+      {
+        if (m_valueLB > val)
+        {
+          builder.AddSegment(ind, leftPos, bottomPos);
+          builder.AddSegment(ind, rightPos, topPos);
+        }
+        else
+        {
+          builder.AddSegment(ind, topPos, leftPos);
+          builder.AddSegment(ind, bottomPos, rightPos);
+        }
+      }
+    }
+  }
+
+  ms::LatLon InterpolatePoint(Square::Rib rib, ValueType alt)
+  {
+    double alt1;
+    double alt2;
+    double lat;
+    double lon;
+
+    switch (rib)
+    {
+    case Rib::Left:
+      alt1 = static_cast<double>(m_valueLB);
+      alt2 = static_cast<double>(m_valueLT);
+      lon = m_left;
+      break;
+    case Rib::Right:
+      alt1 = static_cast<double>(m_valueRB);
+      alt2 = static_cast<double>(m_valueRT);
+      lon = m_right;
+      break;
+    case Rib::Top:
+      alt1 = static_cast<double>(m_valueLT);
+      alt2 = static_cast<double>(m_valueRT);
+      lat = m_top;
+      break;
+    case Rib::Bottom:
+      alt1 = static_cast<double>(m_valueLB);
+      alt2 = static_cast<double>(m_valueRB);
+      lat = m_bottom;
+      break;
+    default:
+      UNREACHABLE();
+    }
+
+    CHECK_NOT_EQUAL(alt, alt2, ());
+    double const coeff = (alt1 - alt) / (alt - alt2);
+
+    switch (rib)
+    {
+    case Rib::Left:
+    case Rib::Right:
+      lat = (m_bottom + m_top * coeff) / (1 + coeff);
+      break;
+    case Rib::Bottom:
+    case Rib::Top:
+      lon = (m_left + m_right * coeff) / (1 + coeff);
+      break;
+    default:
+      UNREACHABLE();
+    }
+
+    return {lat, lon};
+  }
+
+  ValueType m_minValue;
+  ValueType m_valueStep;
+
+  double m_left;
+  double m_right;
+  double m_bottom;
+  double m_top;
+
+  ValueType m_valueLB;
+  ValueType m_valueLT;
+  ValueType m_valueRT;
+  ValueType m_valueRB;
+};
+}  // topography_generator
diff --git a/topography_generator/tile_filter.hpp b/topography_generator/tile_filter.hpp
new file mode 100644
index 0000000000..c63a4a1365
--- /dev/null
+++ b/topography_generator/tile_filter.hpp
@@ -0,0 +1,101 @@
+#pragma once
+
+#include "topography_generator/filters_impl.hpp"
+
+namespace topography_generator
+{
+template <typename ValueType>
+class Filter
+{
+public:
+  virtual ~Filter() = default;
+
+  virtual size_t GetKernelRadius() const = 0;
+  virtual void Process(size_t tileSize, size_t tileOffset,
+               std::vector<ValueType> const & srcValues,
+               std::vector<ValueType> & dstValues) const = 0;
+};
+
+template <typename ValueType>
+class MedianFilter : public Filter<ValueType>
+{
+public:
+  explicit MedianFilter(size_t kernelRadius)
+    : m_kernelRadius(kernelRadius)
+  {}
+
+  size_t GetKernelRadius() const override { return m_kernelRadius; }
+
+  void Process(size_t tileSize, size_t tileOffset,
+               std::vector<ValueType> const & srcValues,
+               std::vector<ValueType> & dstValues) const override
+  {
+    ProcessMedian(m_kernelRadius, tileSize, tileOffset, srcValues, dstValues);
+  }
+
+private:
+  size_t m_kernelRadius;
+};
+
+template <typename ValueType>
+class GaussianFilter : public Filter<ValueType>
+{
+public:
+  GaussianFilter(double standardDeviation, double radiusFactor)
+    : m_kernel(CalculateGaussianLinearKernel(standardDeviation, radiusFactor))
+  {}
+
+  size_t GetKernelRadius() const override { return m_kernel.size() / 2; }
+
+  void Process(size_t tileSize, size_t tileOffset,
+               std::vector<ValueType> const & srcValues,
+               std::vector<ValueType> & dstValues) const override
+  {
+    ProcessWithLinearKernel(m_kernel, tileSize, tileOffset, srcValues, dstValues);
+  }
+
+private:
+  std::vector<double> m_kernel;
+};
+
+template <typename ValueType>
+using FiltersSequence = std::vector<std::unique_ptr<Filter<ValueType>>>;
+
+template <typename ValueType>
+std::vector<ValueType> FilterTile(FiltersSequence<ValueType> const & filters,
+                                  ms::LatLon const & leftBottom,
+                                  size_t stepsInDegree, size_t tileSize,
+                                  ValuesProvider<ValueType> & valuesProvider)
+{
+  size_t combinedOffset = 0;
+  for (auto const & filter : filters)
+    combinedOffset += filter->GetKernelRadius();
+
+  std::vector<ValueType> extTileValues;
+  GetExtendedTile(leftBottom, stepsInDegree, tileSize, combinedOffset, valuesProvider, extTileValues);
+
+  std::vector<ValueType> extTileValues2(extTileValues.size());
+
+  size_t const extTileSize = tileSize + 2 * combinedOffset;
+  CHECK_EQUAL(extTileSize * extTileSize, extTileValues.size(), ());
+
+  size_t currentOffset = 0;
+  for (auto const & filter : filters)
+  {
+    currentOffset += filter->GetKernelRadius();
+    filter->Process(extTileSize, currentOffset, extTileValues, extTileValues2);
+    extTileValues.swap(extTileValues2);
+  }
+
+  std::vector<ValueType> result(tileSize * tileSize);
+  for (size_t i = combinedOffset; i < extTileSize - combinedOffset; ++i)
+  {
+    for (size_t j = combinedOffset; j < extTileSize - combinedOffset; ++j)
+    {
+      size_t const dstIndex = (i - combinedOffset) * tileSize + j - combinedOffset;
+      result[dstIndex] = extTileValues[i * extTileSize + j];
+    }
+  }
+  return result;
+}
+}  // namespace topography_generator
diff --git a/topography_generator/utils/contours.hpp b/topography_generator/utils/contours.hpp
new file mode 100644
index 0000000000..64c3d31a86
--- /dev/null
+++ b/topography_generator/utils/contours.hpp
@@ -0,0 +1,83 @@
+#pragma once
+
+#include "generator/feature_helpers.hpp"
+
+#include "geometry/point2d.hpp"
+#include "geometry/region2d.hpp"
+
+#include <vector>
+#include <unordered_map>
+
+namespace topography_generator
+{
+using Contour = std::vector<m2::PointD>;
+
+template <typename ValueType>
+struct Contours
+{
+  std::unordered_map<ValueType, std::vector<Contour>> m_contours;
+  ValueType m_minValue;
+  ValueType m_maxValue;
+  ValueType m_valueStep;
+  size_t m_invalidValuesCount = 0; // for debug purpose only.
+};
+
+template <typename ValueType>
+void CropContours(m2::RectD & rect, std::vector<m2::RegionD> & regions, size_t maxLength,
+                  Contours<ValueType> & contours)
+{
+  contours.m_minValue = std::numeric_limits<ValueType>::max();
+  contours.m_maxValue = std::numeric_limits<ValueType>::min();
+
+  auto it = contours.m_contours.begin();
+  while (it != contours.m_contours.end())
+  {
+    std::vector<Contour> levelCroppedContours;
+    for (auto const & contour : it->second)
+    {
+      Contour cropped;
+      cropped.reserve(contour.size());
+      for (auto const & pt : contour)
+      {
+        cropped.push_back(pt);
+        if (!rect.IsPointInside(pt) || !RegionsContain(regions, pt) || cropped.size() == maxLength)
+        {
+          if (cropped.size() > 1)
+            levelCroppedContours.emplace_back(std::move(cropped));
+          cropped = {};
+        }
+      }
+      if (cropped.size() > 1)
+        levelCroppedContours.emplace_back(std::move(cropped));
+    }
+    it->second.swap(levelCroppedContours);
+
+    if (!it->second.empty())
+    {
+      contours.m_minValue = std::min(it->first, contours.m_minValue);
+      contours.m_maxValue = std::max(it->first, contours.m_maxValue);
+      ++it;
+    }
+    else
+    {
+      it = contours.m_contours.erase(it);
+    }
+  }
+}
+
+template <typename ValueType>
+void SimplifyContours(int simplificationZoom, Contours<ValueType> & contours)
+{
+  for (auto & levelContours : contours.m_contours)
+  {
+    for (auto & contour : levelContours.second)
+    {
+      std::vector<m2::PointD> contourSimple;
+      feature::SimplifyPoints(m2::SquaredDistanceFromSegmentToPoint<m2::PointD>(),
+                              simplificationZoom, contour, contourSimple);
+      contour.swap(contourSimple);
+    }
+  }
+}
+}  // namespace topography_generator
+
diff --git a/topography_generator/utils/contours_serdes.hpp b/topography_generator/utils/contours_serdes.hpp
new file mode 100644
index 0000000000..a03948b9e1
--- /dev/null
+++ b/topography_generator/utils/contours_serdes.hpp
@@ -0,0 +1,138 @@
+#pragma once
+
+#include "topography_generator/utils/contours.hpp"
+
+#include "coding/file_writer.hpp"
+#include "coding/file_reader.hpp"
+#include "coding/geometry_coding.hpp"
+
+namespace topography_generator
+{
+template <typename ValueType>
+class SerializerContours
+{
+public:
+  explicit SerializerContours(Contours<ValueType> && contours): m_contours(std::move(contours)) {}
+
+  template <typename Sink>
+  void Serialize(Sink & sink)
+  {
+    WriteToSink(sink, m_contours.m_minValue);
+    WriteToSink(sink, m_contours.m_maxValue);
+    WriteToSink(sink, m_contours.m_valueStep);
+
+    WriteToSink(sink, static_cast<uint32_t>(m_contours.m_contours.size()));
+    for (auto const & levelContours : m_contours.m_contours)
+    {
+      SerializeContours(sink, levelContours.first, levelContours.second);
+    }
+  }
+private:
+  template <typename Sink>
+  void SerializeContours(Sink & sink, ValueType value,
+                         std::vector<topography_generator::Contour> const & contours)
+  {
+    WriteToSink(sink, value);
+    WriteToSink(sink, static_cast<uint32_t>(contours.size()));
+    for (auto const & contour : contours)
+    {
+      SerializeContour(sink, contour);
+    }
+  }
+
+  template <typename Sink>
+  void SerializeContour(Sink & sink, topography_generator::Contour const & contour)
+  {
+    WriteToSink(sink, static_cast<uint32_t>(contour.size()));
+    sink.Write(contour.data(), contour.size() * sizeof(contour[0]));
+    /*
+    serial::GeometryCodingParams codingParams(kFeatureSorterPointCoordBits, 0);
+    codingParams.SetBasePoint(contour[0]);
+    std::vector<m2::PointD> toSave;
+    toSave.insert(contour.begin() + 1, contour.end());
+    serial::SaveInnerPath(toSave, codingParams, sink);
+    */
+  }
+
+  Contours<ValueType> m_contours;
+};
+
+template <typename ValueType>
+class DeserializerContours
+{
+public:
+  template <typename Reader>
+  void Deserialize(Reader & reader, Contours<ValueType> & contours)
+  {
+    NonOwningReaderSource source(reader);
+    contours.m_minValue = ReadPrimitiveFromSource<ValueType>(source);
+    contours.m_maxValue = ReadPrimitiveFromSource<ValueType>(source);
+    contours.m_valueStep = ReadPrimitiveFromSource<ValueType>(source);
+
+    size_t const levelsCount = ReadPrimitiveFromSource<uint32_t>(source);
+    for (size_t i = 0; i < levelsCount; ++i)
+    {
+      ValueType levelValue;
+      std::vector<topography_generator::Contour> levelContours;
+      DeserializeContours(source, levelValue, levelContours);
+      contours.m_contours[levelValue].swap(levelContours);
+    }
+  }
+
+private:
+  void DeserializeContours(NonOwningReaderSource & source, ValueType & value,
+                            std::vector<topography_generator::Contour> & contours)
+  {
+    value = ReadPrimitiveFromSource<ValueType>(source);
+    size_t const contoursCount = ReadPrimitiveFromSource<uint32_t>(source);
+    contours.resize(contoursCount);
+    for (auto & contour : contours)
+      DeserializeContour(source, contour);
+  }
+
+  void DeserializeContour(NonOwningReaderSource & source,
+                          topography_generator::Contour & contour)
+  {
+    size_t const pointsCount = ReadPrimitiveFromSource<uint32_t>(source);
+    contour.resize(pointsCount);
+    source.Read(contour.data(), pointsCount * sizeof(contour[0]));
+  }
+};
+
+template <typename ValueType>
+bool SaveContrours(std::string const & filePath,
+                   Contours<ValueType> && contours)
+{
+  try
+  {
+    FileWriter file(filePath);
+    SerializerContours<ValueType> ser(std::move(contours));
+    ser.Serialize(file);
+  }
+  catch (FileWriter::Exception const & ex)
+  {
+    LOG(LWARNING, ("File writer exception raised:", ex.what(), ", file", filePath));
+    return false;
+  }
+  return true;
+}
+
+template <typename ValueType>
+bool LoadContours(std::string const & filePath, Contours<ValueType> & contours)
+{
+  {
+    try
+    {
+      FileReader file(filePath);
+      DeserializerContours<ValueType> des;
+      des.Deserialize(file, contours);
+    }
+    catch (FileReader::Exception const & ex)
+    {
+      LOG(LWARNING, ("File writer exception raised:", ex.what(), ", file", filePath));
+      return false;
+    }
+    return true;
+  }
+}
+}  // namespace topography_generator
diff --git a/topography_generator/utils/values_provider.hpp b/topography_generator/utils/values_provider.hpp
new file mode 100644
index 0000000000..0d3664ae45
--- /dev/null
+++ b/topography_generator/utils/values_provider.hpp
@@ -0,0 +1,14 @@
+#pragma once
+
+#include "geometry/latlon.hpp"
+
+namespace topography_generator
+{
+template <typename ValueType>
+class ValuesProvider
+{
+public:
+  virtual ValueType GetValue(ms::LatLon const & pos) = 0;
+  virtual ValueType GetInvalidValue() const = 0;
+};
+}  // namespace topography_generator