Tomasz Bartuś

Atrybut funkcyjny w GeoMedia obliczający liczbę punktów o danej wartości atrybutu

Tomasz Bartuś

OPROGRAMOWANIE:

Geomedia Professional: 06.01.11.13

PROBLEM:

Siatka kwadratowa o boku 500 m.

Wewnątrz oczek siatki znajdują się obiekty punktowe o różnej wartości atrybutu ranga. Dla każdego oczka siatki mamy wykonać proste obliczenie arytmetyczne dane wzorem:

gdzie:

m - liczba wszystkich jednostek w danym oczku siatki,

np: w przypadku przedstawionym na Fig. 1:

1) Oczko lewe górne (I):
liczba wszystkich jednostek w oczku siatki I: 6;
liczba elementów zielonych (A): 2;
liczba elementów czerwonych (B): 1;
liczba elementów niebieskich (C): 3;

H_I=-Σ(q_I*ln q_I)= -[(2/6*ln2/6)+(1/6*ln1/6)+(3/6*ln3/6)] 2) Oczko prawe górne (II):
liczba wszystkich jednostek w oczku siatki II: 0;
liczba elementów zielonych (A): 0;
liczba elementów czerwonych (B): 0;
liczba elementów niebieskich (C): 0;

H_II=-Σ(q_II*ln q_II)=0 3) Oczko lewe dolne (III):
liczba wszystkich jednostek w oczku siatki III: 1;
liczba elementów zielonych (A): 1;
liczba elementów czerwonych (B): 0;
liczba elementów niebieskich (C): 0;

H_III=-Σ(q_III*ln q_III)= -[(1/1*ln1/1)]=0 4) Oczko prawe dolne (IV):
liczba wszystkich jednostek w oczku siatki IV: 0;
liczba elementów zielonych (A): 0;
liczba elementów czerwonych (B): 0;
liczba elementów niebieskich (C): 0;

H_IV=-Σ(q_IV*ln q_IV)=0

ROZWIĄZANIE:

[dziękuję za olbrzymią pomoc P. Leszkowi Nowickiemu]:

1. W oknie dialogowym Feature Class, otwieramy kolejno klasę punktów (np: geologia_punktowe) i klasę siatki (np: siatka) i zmieniamy mazwy ich kluczy podstawowych na intuicyjne, kolejno: ID_geologia_punkty oraz ID_siatka. Krok ten nie jest konieczny ale ułatwi nam późniejsze prace.
2. Obliczamy iloczyn przestrzenny Spartial intersection klasy punktów (geologia_punktowe) i klasy siatki (siatka) (zapytanie wynikowe: geol_pkt_siatka). Każdy z punktów otrzymuje w ten sposób dodatkowy atrybut opisujący jego położenie względem oczek siatki (ID_siatka).
3. Za pomoca polecenia Aggregation, dla każdego okna siatki, dokonujemy obliczeń liczby wszystkich elementów. Załóżmy, że mamy 8 typów obiektów punktowych (rang) (np: kamieniołomy, piaskownie, jaskinie, ostańce denudacyjne itd.). Sumujemy najpierw obiekty danej rangi, a później sumy wszystkich rang. Należy zauważyć, że w podanym poniżej wzorze wyeliminowano obiekty punktowe o randze "4" i nie biorą one udziału w obliczeniach. Aggregate to summary features in: siatka; From detail features in: geol_pkt_siatka; Spartial Aggregation: touch; Attribute Aggregation: Id_siatka/Id_siatka; Output: NEW. Tworzymy atrybut funkcyjny dla każdego okna siatki obliczający liczbę wszystich punktów, tzn. m (L_j):
   
   SUM(SUM(IF(Detail.ranga=1; 1; 0)); SUM(IF(Detail.ranga=2; 1; 0)); SUM(IF(Detail.ranga=3; 1; 0)); SUM(IF(Detail.ranga=5; 1; 0)); SUM(IF(Detail.ranga=6; 1; 0)); SUM(IF(Detail.ranga=7; 1; 0)); SUM(IF(Detail.ranga=8; 1; 0)))
   
   Utworzony atrybut funkcyjny zapisujemy pod nazwą: SUMAelementów_Lj.
4. Przechodzimy do kolejnego kroku, t.j. obliczenia dla każdego oczka siatki, liczby elementów danego typu (rangi). Jeżeli nie wychodziliśmy z okna Aggregation, wystarczy w zakładce: Output; obok okna Output functional attributes, kliknąć new i wprowadzić formułę: SUM(IF(Detail.ranga=1; 1; 0))
   zapisując ją pod nową nazwą atrybutu funkcyjnego np: SUMA1. To samo robimy dla kolejnych rang punktów:
   
   SUM(IF(Detail.ranga=2; 1; 0)) -> nazwa: SUMA2
   SUM(IF(Detail.ranga=3; 1; 0)) -> nazwa: SUMA3
   SUM(IF(Detail.ranga=5; 1; 0)) -> nazwa: SUMA5
   SUM(IF(Detail.ranga=6; 1; 0)) -> nazwa: SUMA6
   SUM(IF(Detail.ranga=7; 1; 0)) -> nazwa: SUMA7
   SUM(IF(Detail.ranga=8; 1; 0)) -> nazwa: SUMA8
   
   Utworzone atrybuty funkcyjne, co prawda nie rozwiązują zasadniczego problemu ale staną się pomocne przy testowaniu poprawności głównej formuły.
5. Przechodzimy do najważniejszej cześci, t.j. opracowania atrybutu funkcyjnego dla każdego oczka siatki obliczającego wartość parametru H_jed. Możemy przyjąć, że formuła będzie miała następującą postać ogólną:
   
   -SUM(
   obliczenia dla rangi 1;
   obliczenia dla rangi 2;
   obliczenia dla rangi 3;
   obliczenia dla rangi 5;
   obliczenia dla rangi 6;
   obliczenia dla rangi 7;
   obliczenia dla rangi 8;
   )
   
   Tworzymy obliczenia dla rangi 1. W obliczeniach wykorzystujemy działanie utworzonego wcześniej atrybutu: SUMAelementów_Lj
   
   (CAST(SUM(IF(Detail.ranga=1; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=1; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=1; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
   
   i dalej kolejno:
   
   • (CAST(SUM(IF(Detail.ranga=2; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=2; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=2; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
   • (CAST(SUM(IF(Detail.ranga=3; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=3; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=3; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
   • (CAST(SUM(IF(Detail.ranga=5; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=5; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=5; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
   • (CAST(SUM(IF(Detail.ranga=6; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=6; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=6; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
   • (CAST(SUM(IF(Detail.ranga=7; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=7; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=7; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
   • (CAST(SUM(IF(Detail.ranga=8; 1; 0))/Output.SUMAelementów_Lj; double))*(IF(SUM(IF(Detail.ranga=8; 1; 0))>0; (LN(CAST(SUM(IF(Detail.ranga=8; 1; 0))/Output.SUMAelementów_Lj; Double))); 0))
6. Łącząc utworzone elementy w sumę, tworzymy atrybut funkcyjny obliczający H_jed:
   
   -SUM(CAST(SUM(IF(Detail.ranga=1; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=1; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=1; 1; 0))/Output.SUMAelementów_Lj; Double)); 0); CAST(SUM(IF(Detail.ranga=2; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=2; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=2; 1; 0))/Output.SUMAelementów_Lj; Double)); 0); CAST(SUM(IF(Detail.ranga=3; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=3; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=3; 1; 0))/Output.SUMAelementów_Lj; Double)); 0); CAST(SUM(IF(Detail.ranga=5; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=5; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=5; 1; 0))/Output.SUMAelementów_Lj; Double)); 0); CAST(SUM(IF(Detail.ranga=6; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=6; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=6; 1; 0))/Output.SUMAelementów_Lj; Double)); 0); CAST(SUM(IF(Detail.ranga=7; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=7; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=7; 1; 0))/Output.SUMAelementów_Lj; Double)); 0); CAST(SUM(IF(Detail.ranga=8; 1; 0))/Output.SUMAelementów_Lj; Double)*IF(SUM(IF(Detail.ranga=8; 1; 0))>0; LN(CAST(SUM(IF(Detail.ranga=8; 1; 0))/Output.SUMAelementów_Lj; Double)); 0)) Atrybut nazywamy: Hjed.
7. Na koniec, choć to oczywiście nie jest niekonieczne, dla odróżnienia sytuacji, w której w oczku siatki brak jakichkolwiek elementów punktowych (H_jed=0), od sytuacji, w której mamy dowolną ilość punktów tej samej rangi (H_jed=0), możemy dołożyć warunek:
   
   IF(Output.SUMAelementów_Lj=0; 99; formuła atrybutu funkcyjnego Hjed)
   
   Teraz w oczkach siatki, w których brak jakichkolwiek elementów, w wyniku działania atrybutu, pojawi się wartość "99".