Documentação

osmc:Metodologia/Algoritmo SQL/Issues: mudanças entre as edições

osmc Discussão

osmc:Metodologia/Algoritmo SQL/Issues (ver código-fonte)

Edição das 23h47min de 25 de junho de 2024

13 228 bytes adicionados ,  25 junho
m
→‎Issue 05 - Cálculo de dígitos base32: typo
(Criou página com 'Issues e discussões conceituais, antes de se partir para a implementação. == Issue 01 - Uso ótimo do hInt64 == Nesta issue buscamos garantir melhor ocupação nos 57 disponíveis de geocódigo hInt64. Os blocos hierárquicos de bits seriam os seguintes, concatenados: : <code>$country$L0$grid</code> * O bloco <code>$grid</code> não pode ser otimizado. Possui como principal requsito a garantia de um número par de b...')
 
 
(16 revisões intermediárias pelo mesmo usuário não estão sendo mostradas)
Linha 1: Linha 1:
Issues e discussões conceituais, antes de se partir para a implementação.
Issues e discussões conceituais, antes de se partir para a implementação. Ver também [[Discussão:DNGS/Decisões_soberanas#Issues]], relativas a decisões humanas, sem algoritmo previsto.


== Issue 01 - Uso ótimo do hInt64 ==
== Issue 01 - Uso ótimo do hInt64 ==
Linha 17: Linha 17:
=== Bloco L0 ===
=== Bloco L0 ===


Se a cobertura do país for exatamente uma célula quadrada, pode ser ausente, mas em geral sua geometria precisa ser otimizada, para mais de 70% da área de L0 ser ocupada pelo polígono do território.
Se a cobertura do país for exatamente uma célula quadrada, pode ser ausente, mas em geral sua geometria precisa ser otimizada, para mais de 45% da área de L0 ser ocupada pelo polígono do território. Sem essa exigência de ocupação as aplicações de abrangência  global perdem performance. As interseções entre países precisam ser minimizadas.


O máximo de células é 16, o mínimo 2 (em países com fronteiras instáveis o mais seguro é forçar reserva). Em países com ilhas e células de cobertura com baixo percentual de uso (e possibilidade de representação com uma só célula-filha) pode-se adotar mais que 16 células de cobertura, emulando um nível L0 com mosaico de pedaços, como no caso do Brail. Vamos discutir os três casos a seguir. <br/> A questão principal é se vai permitir ou não quantidade ímpar de células.
<syntaxhighlight lang="sql">
WITH m as (
  select pais, case when count(*)>16 then 16 else count(*) end*max(area_km2) as mx
  from report.v001_osmc_coverage_l0_list group by 1
)
SELECT v.pais, count(*) n, sum(v.area_km2) as km2, round(100.0*sum(v.area_km2/m.mx))::text||'%' perc_coberto
FROM report.v001_osmc_coverage_l0_list v INNER JOIN m ON m.pais=v.pais
GROUP BY 1;
</syntaxhighlight>
<pre>
pais| n  |  km2  | perc_coberto
BR | 18 | 8729499 | 50%
CM | 14 |  466041 | 48%
CO | 16 | 1886770 | 44%
</pre>


No caso do Brasil a codificação L0 ...
O máximo de células é 16, o mínimo 2 (em países com fronteiras instáveis o mais seguro é forçar reserva). Em países com ilhas e células de cobertura com baixo percentual de uso (e possibilidade de representação com uma só célula-filha) pode-se adotar mais que 16 células de cobertura, emulando um nível L0 com mosaico de pedaços, como no caso do Brail. Vamos discutir os três casos a seguir. <br /> A questão principal é se vai permitir ou não quantidade ímpar de células.
 
No caso do Brasil a codificação L0 requer 2 conjuntos de ilhas, depois das 16 coberturas territoriais. Basta 1 bit para codificar se é uma "falsa célula" ou não. O extremo sul do Brasil, que recebeu o rótulo <code>f</code>, teve sua área subdivida em mosaico para comportar as "falsas células" <code>10=g</code> e  <code>11=h</code>. Na prática são todas  <code>f</code>, e o geocódigo só "assume sua representação concreta" a partir do segundo dígito: <code>fT</code> (e seu litoral <code>fY</code>) para o território extremo sul, <code>fN</code> para as ilhas oceânicas sob jurisdição do estado de Espirito Santo, e <code>fP</code> para o arquipelogo de São Pedro e São Paulo.
 
[[Arquivo:BR-L0coverSci-f.png|miniaturadaimagem|420px|Goecódigo científico do Brasil. Geometria concreta das células-filhas de "f", partida em 8. Apenas 4 partes em uso, restante são reservas. A grade existe para as células <code>fT</code> e <code>fY</code> no mesmo local, e para as coberturas oceânicas <code>fP</code> e <code>fN</code>.]]
 
{|
|-
|<pre>
pais | b16h | area_km2
------+------+----------
BR  | 00  |  275241
BR  | 01  |  757637
BR  | 02  |  645014
BR  | 03  |  132780
BR  | 04  |  619969
BR  | 05  |  1075008
BR  | 06  |  1099512
BR  | 07  |  1054743
BR  | 08  |  103225
</pre>
|<pre>
pais | b16h | area_km2
------+------+----------
BR  | 09  |  333725
BR  | 0a  |  1099383
BR  | 0b  |  723005
BR  | 0c  |    19403
BR  | 0d  |  714206
BR  | 0e  |    54071
BR  | 0f  |    17440
BR  | 10  |    1519
BR  | 11  |    3618
</pre>
|}
 
== Issue 02 - Função estimadora de cobertura interior e quadrado envolvente ==
 
Para fazer estimativas sobre células de km2. A função f(x) que retorna potências ''side'' e ''next'':
 
* ''side'': estimativa para conferir se a área de uma célula recortada (ou seja célula de cobertura) é ou não uma célula interior (sem recortes), apenas pelo valor de sua área.
* ''next'': para chutar o quadrado envolvente inicial de um país (pode dispensar cobertura), ou de um município.
 
Abaixo a função com exemplo da área do Brasil.
 
<pre>
-- A = 8510417.771 km2:
with t as (select sqrt(8510417.771 *1000^2) as x)
select p as side_pow_exp, 2^p as side_pow_val,
      round( (100.0*(2^p-x)/2^p)::numeric , 2)::real  side_diff_perc,
      p+1 as next_pow_exp, 2^(p+1) as next_pow_val,
      round( (100.0*(2^(p+1)-x)/2^(p+1))::numeric , 2)::real  next_diff_perc
from (select round(log(x)/log(2)) p,x  from t) t2;
-- side_pow_exp | side_pow_val | side_diff_perc | next_pow_exp | next_pow_val | next_diff_perc
-- -------------+--------------+----------------+--------------+--------------+----------------
--          21 |      2097152 |        -39.11 |          22 |      4194304 |          30.45
-- A = 757637 km2:
side_pow_exp | side_pow_val | side_diff_perc |
--------------+--------------+----------------|
          20 |      1048576 |          16.99 |
-- A = 1099512 km2:
          20 |      1048576 |              0 |
</pre>
No caso da área do Brasil (de 8510417.771 km2) temos o next como estimativa de quadrado envolvente, com lados de 2^22 km:
# quadrado com 2^22 metros de lado contém o Brasil.
# Primeiro dígito base4 é a cobertura de 4 quadrados com 2^21 metros cada
# Primeiro dígito base16 é a cobertura de 16 quadrados com 2^20 metros de lado cada.
 
Nas áreas de células conferimos qual tem ''side_diff_perc=0'', que foi apenas o caso de A = 1099512 km2.
 
Exemplos, onde seriam esperados expoentes pares,  da ordem de 20 para cobertura do Brasil e 18 para Camarões.
 
{| class="wikitable"
!        pais          !! is_contained !!        cbits        !! b16 !! area_km2 !!  exp e diff%
|-
| BR = 76 = 0001001100  || f            || 0001001100.00000110 || 06  ||  1099512 ||  20 com 0%
|-
| BR = 76 = 0001001100  || t            || 0001001100.00000000 || 00  ||  275241 ||  19 com -0.07%
|-
| BR = 76 = 0001001100  || t            || 0001001100.00000001 || 01  ||  757637 ||  20 com 16.99%
|-
|...
|-
| CM = 120 = 0001111000 || f            || 0001111000.1001    || 9    ||    68719 || 18 com 0%   
|-
| CM = 120 = 0001111000 || t            || 0001111000.0001    || 1    ||    2469 || 16 com 24.18%
|-
|...
|}
 
A interseção da geometria da célula com o polígono do país tem a sua área indicada pela coluna ''area_km2''. A coluna ''is_contained'' é um flag, verdadeiro quando a célula tem interseção com as bordas do país. Analisando as áreas do Brasil, o ''flag'' foi falso apenas para os <math>1099512 \ km^2</math> de área, cuja raíz quadada é o tamanho do lado do quadrado, e uma das [[osmc:Metodologia#Projeção_e_cobertura_nacionais|potências de 2 utilizadas]], <math>\sqrt{1099512\ km^2} = 1048576 \ m = 2^{20} \ m = 1048,58 \ km</math>. Analogamente em Camarões (CM),  <math>\sqrt{68719\ km^2} = 262143.1 \ m \approx 262144 \ m = 2^{18} \ m = 262.144 \ km</math>.
 
== Issue 03 - Redução do excesso de pontos por efeito-varal  ==
Conforme [https://gis.stackexchange.com/q/444441/7505 "How to reduce the clothesline-effect in ST_Transform?"], pode-se minimizar o efeito-varal. A ideia aqui nesta issue é usasar simplificação de curva, para controle e eventual complementação do algoritmo.
 
<syntaxhighlight lang="sql">
select hlevel, count(*) n_cells, round(sqrt(max(st_area(geom)))/1000.0) side_km,
      round(avg(ST_NPoints(geom))) geom_npts,
      round(avg(ST_CharactDiam(geom))) as geom_chdiam, -- relative ref.
      round(avg(ST_NPoints(geom4326))) g4326_npts,
      round(avg(ST_NPoints(st_simplify(geom4326,0.00000005)))) g4326simpl_npts
      -- 0.00000005 is SRID-metric absolute ref
from grid_cm.all_levels
group by 1 order by 1;
</syntaxhighlight>
 
<pre>
hlevel | n_cells | side_km | npts_geom | npts_geom | npts_simpl
--------+---------+---------+-----------+-----------+------------
      0 |      14 |    262 |        5 |      749 |        293
    0.5 |      28 |    185 |        5 |      769 |        280
      1 |      56 |    131 |        5 |      749 |        152
    1.5 |    112 |      93 |        5 |      769 |        141
      2 |    224 |      66 |        5 |      749 |        77
    2.5 |    448 |      46 |        5 |      769 |        71
      3 |    896 |      33 |        5 |      749 |        39
    3.5 |    1792 |      23 |        5 |      769 |        36
      4 |    3584 |      16 |        5 |      749 |        20
    4.5 |    7168 |      12 |        5 |      769 |        19
      5 |  14336 |      8 |        5 |      749 |        11
    5.5 |  28672 |      6 |        5 |      769 |        11
</pre>
 
Conclusão: refazer a nossa função de ''ST_Transform_resilient()'' para os casos de cache, talvez não de alta performance
# no final incluir ST_SimplifyPreserveTopology, portanto mais esse parametro ''simplify_tolerance'' precisa ser incluso. Se o uso se restringe a células, pode ser ST_Simplify, mais rápida. Decidir se existe como criar uma tolerância relativa ao perímetro, abarcando a intenção da densidade. Provavelmente um fator da densidade*perimetro resolva. O foco são as células da grade.
# com base no parâmetro métrico  ''simplify_tolerance'' decidir se não faz nada, visto que depois do simplify voltaria ao número de pontos original. É de fato o esperado para a grade de hlevel de 6 em diante no exemplo. Foco nas células.
 
Na resposta da questão quebrar em duas,  "general" depois "grids" (exemplificar com níveis inteiros pulando os meios).  Nas grades supor que só um eixo precisa de pontos ... Ver eles no QGIS ou pelo menos nos dados de 11 pontos. O segmentize aplicado a uma só direção é complexo, não vale a pena.
 
<syntaxhighlight lang="sql">
DROP FUNCTION IF EXISTS ST_Transform_Resilient
;
CREATE FUNCTION ST_Transform_Resilient(
  g    geometry,    -- the input geometry
  srid integer,      -- target SRID, to transform g
  size_fraction float DEFAULT 0.05,  -- 1/density. Density of points per charactDiam (or negative for absolute fraction).
  tolerance    float DEFAULT 0      -- E.g. on srid=4326 use 0.00000005.
                                    -- ZERO=0.000000001. Smallest=0.000000005. Start cell with 0.00000002.
) RETURNS      geometry
  language SQL IMMUTABLE
AS $f$
SELECT CASE
      WHEN COALESCE(size_fraction,0.0)>0.0 AND COALESCE(tolerance,0)>0 THEN
          ST_SimplifyPreserveTopology(geom,tolerance) -- ST_Simplify enough for grid cells
      ELSE geom
      END
FROM (
  SELECT CASE
    WHEN size>0.0 THEN  ST_Transform(  ST_Segmentize(g,size)  , srid  )
    ELSE  ST_Transform(g,srid)
    END geom, size
  FROM (
    SELECT CASE
        WHEN size_fraction IS NULL THEN 0.0
        WHEN size_fraction<0 THEN -size_fraction
        ELSE ST_CharactDiam(g) * size_fraction
        END
  ) t1(size)
) t2
$f$;
COMMENT ON FUNCTION ST_Transform_Resilient IS
  'Avoid clothesline-effect. See problem/solution discussed at https://gis.stackexchange.com/q/444441/7505'
;
</syntaxhighlight>
 
Migrar função para o git oficial!
 
== Issue 04 - Lib de grades  ==
Desenvolvendo em [[osmc:Metodologia/Algoritmo SQL/Issue04]].
 
== Issue 05 - Cálculo de dígitos base32 ==
 
Por imposição do [[Discrete National Grid Systems/pt|padrão DNGS]] temos <math>S_{L}=2^{Lmax-L}</math>. No caso do Brasil ''Lmax''=20, no caso de Camarões ''Lmax=18''. Já havíamos notado isso ao apresentar a [[osmc:Metodologia/Algoritmo_SQL/Lib#Core|Lib Core]].
 
No caso de geocódigo logístico, que requer cobertura municipal, supondo que o nível da grade adotada para cobertura seja ''Lcover'', podemos obter a aproximação ao metro por um múltiplo de 2.5, que é a quantidade de níveis que se sobe a cada 5 bits do dígito da base32. Partindo de  <math>1=2^{Lmax - Lcover - N \cdot 2.5}</math> obtemos por <math>log_2(1)=0</math> a fórmula <math>N= (Lmax-Lcover)/2.5</math>, onde ''N'' é o número de dígitos depois do primeiro para se chegar ao metro.  Num caso típico a cobertura municipal no Brasil (''Lmax''=20) é de nível ''L7.5'', resultando em ''N''=5, ou seja, com 6 dígitos teremos o metro. Para uma cobertura de resolução menor, ''L7'', o valor já não é exato, ''N''=5.2. O número de dígitos para chegar no metro portanto será <math>Ndig1m = 1+\lceil (Lmax-Lcover)/2.5 \rceil</math>.
 
<syntaxhighlight lang="sql" style="font-size: 80%;">
CREATE FUNCTION grid_br.logistic_Ndig1m(intlevel_cover int) RETURNS real[] AS $f$
    SELECT array[  Lcover, ceil(Ndig1m), ceil(Ndig1m)-Ndig1m,
          2^(Lmax-Lcover-ceil(Ndig1m)*2.5),
          CASE WHEN ceil(Ndig1m)-Ndig1m=0 THEN NULL ELSE 2^(Lmax-Lcover-floor(Ndig1m)*2.5) END  ]
    FROM (
      SELECT *, (Lmax-Lcover)/2.5 as Ndig1m
      FROM (SELECT 20 as Lmax, intlevel_cover/10.0 as Lcover) t1 
    ) t2
$f$ language SQL IMMUTABLE;
 
-- gerando a tabela do Brasil:
select x[1] as level_cover, x[2] as Ndigits_for_1m, round(x[4],2) as side_size,
      round(x[5],2) as size_before
from (
  select grid_br.logistic_Ndig1m(intLevel) as x
  from generate_series(40,90,5) t1(intLevel)
) t2;
</syntaxhighlight>
 
A escolha de nivel ''Lcover'' de cada município no caso do Brasil vai ter o seguinte perfil dado pela função <code>logistic_Ndig1m()</code>. Como o metro é a maior resolução significativa, é interesante não permitir resoluções finais inferiores a meio metro: por exemplo quando a cobertura ótima ''Lcover'' for 4 ou 4.5, usar no máximo 6 dígitos ao invés de 7, já que 2 metros e 1.4 metros são resoluções razoáveis para endereços. Idem ''Lcover''6.5 e 7 etc.
 
{| class="wikitable"
|+ Número de dígitos base32 para se chegar no metro, conforme ''Lcover'' do município
|-
|''Lcover'' || Ndigits for 1m || ''side_size'' || ''size_before''
|-
|        4 ||              7 ||      0.35 ||        2.00
|-
|      4.5 ||              7 ||      0.25 ||        1.41
|-
|        5 ||              6 ||      1.00 ||           
|-
|      5.5 ||              6 ||      0.71 ||        4.00
|-
|        6 ||              6 ||      0.50 ||        2.83
|-
|      6.5 ||              6 ||      0.35 ||        2.00
|-
|        7 ||              6 ||      0.25 ||        1.41
|-
|      7.5 ||              5 ||      1.00 ||           
|-
|        8 ||              5 ||      0.71 ||        4.00
|-
|      8.5 ||              5 ||      0.50 ||        2.83
|-
|        9 ||              5 ||      0.35 ||        2.00
|}
Peter
2 384

edições

Disponível em “https://wiki.addressforall.org/doc/Especial:MobileDiff/3879...4099”