Re: [GRASS-user] compare a DCELL and FCELL question
Markus Metz <[email protected]> Sun, 24 Jan 2021 18:53:44 +0100
| Newsgroups | gmane.comp.gis.grass.devel,gmane.comp.gis.grass.user |
|---|---|
| Message-ID | <CAG+h=FGjb0gEbZ+McR6Jde6RC0QzDygjttVEWS+FAswT2M0Fxg@mail.gmail.com> |
--===============6778851474218773543== Content-Type: multipart/alternative; boundary="000000000000ea919305b9a91aab" --000000000000ea919305b9a91aab Content-Type: text/plain; charset="UTF-8" Content-Transfer-Encoding: quoted-printable A note about the usage of GRASS_EPSILON =3D 1.0e-15: The range of single precision floating point is about 1.4012984643e=E2=88= =9245 to 3.4028234664e38 That means for large numbers, using GRASS_EPSILON will also detect floating point precision limits, not only meaningful differences. For small numbers, GRASS_EPSILON will not detect meaningful differences. GRASS_EPSILON could instead be modified with something like "max(abs(map_A), abs(map_B)) * 1.0e-15" to test for meaningful differences. Markus M On Sun, Jan 24, 2021 at 5:32 PM Markus Metz <[email protected]> wrote: > > > > On Sun, Jan 24, 2021 at 5:06 PM ming han <[email protected]> wrote: > > > > Hi Everyone > > > > Many thanks for your help. Is if(fabs(map_A - map_B) <=3D 1.0e-15, .= .. ) approach may increase the runtime compare to '=3D=3D' way? > > The formulas are very simple, I don't think that differences in runtime can be reliably measured. (De-)compression of the data and the operating system's file cache have a much stronger influence on the runtime. > > Markus M > > > > Thanks > > Ming > > > > Markus Metz <[email protected]> =E4=BA=8E2021=E5=B9=B41=E6= =9C=8824=E6=97=A5=E5=91=A8=E6=97=A5 =E4=B8=8A=E5=8D=8810:57=E5=86=99=E9=81= =93=EF=BC=9A > >> > >> Trying to answer the original question: with a DCELL map "cat1_acc_riv" and a FCELL map "cat1_minacc", why is "float(cat1_acc_riv) =3D=3D float(cat1_minacc)" not equal to "int(cat1_acc_riv) =3D=3D int(cat1_= minacc)" ? > >> > >> int truncates to integer while float converts to single precision floating point. E.g. with cat1_acc_riv =3D 1.1 and cat1_minacc =3D 1.9, "float(cat1_acc_riv) =3D=3D float(cat1_minacc)" becomes "1.1 =3D=3D 1.9" wh= ereas "int(cat1_acc_riv) =3D=3D int(cat1_minacc)" becomes "1 =3D=3D 1", thus the = results are different. > >> > >> Another reason for possible differences is that float can only represent max 7 decimal digits. E.g. float(194320567) becomes 194320560 but int(194320567) preserves the value 194320567. > >> > >> Thus the safest is to cast everything to the type with the highest precision. In this case with FCELL and DCELL, use "double(cat1_acc_riv) =3D= =3D double(cat1_minacc)" or even better the suggestion of Markus N. > >> > >> Markus M > >> > >> > >> On Sun, Jan 24, 2021 at 3:51 PM ming han <[email protected]> wrote: > >> > > >> > Hi Markus and Micha > >> > > >> > I am just trying to find grids have the same values in these two rasters, I will try the threshold approach. > >> > > >> > Thanks > >> > Ming > >> > > >> > Markus Neteler <[email protected]> =E4=BA=8E2021=E5=B9=B41=E6=9C=882= 4=E6=97=A5=E5=91=A8=E6=97=A5 =E4=B8=8A=E5=8D=886:58=E5=86=99=E9=81=93=EF=BC= =9A > >> >> > >> >> Hi Ming, > >> >> > >> >> On Sun, Jan 24, 2021 at 10:49 AM ming han <[email protected]> wrote: > >> >> > > >> >> > Hi Micha > >> >> > > >> >> > Many thanks for your reply. > >> >> > Here is the command I am using: > >> >> > > >> >> > if(float(cat1_acc_riv) =3D=3D float(cat1_minacc), str_r, nul= l()) > >> >> > > >> >> > The str_r is a CELL raster. the result is different when I change it to: > >> >> > if(int(cat1_acc_riv) =3D=3D int(cat1_minacc), str_r, null(= )) > >> >> > >> >> Note that numerical "equality" is better tested with a threshold test > >> >> against the map pixel difference. > >> >> As the threshold, we use GRASS_EPSILON which is defined as 1.0e-15. > >> >> > >> >> Hence the test needs to be implemented in a different way, i.e. by > >> >> using an epsilon. > >> >> Essentially something like this: > >> >> > >> >> if(fabs(map_A - map_B) <=3D 1.0e-15, ... ) > >> >> > >> >> In your case (untested): > >> >> r.mapcalc diffepsilon =3D if( abs( map_A - map_B) <=3D 1.0e-15, str= _r , null()) > >> >> > >> >> See related discussions here: [1], [2] and elsewhere. > >> >> > >> >> [1] Comment by Glynn: https://trac.osgeo.org/grass/ticket/2854#comment:9 > >> >> [2] Comment by Glynn: > >> >> https://lists.osgeo.org/pipermail/grass-user/2015-October/073200.html > >> >> > >> >> Best, > >> >> Markus > >> > > >> > _______________________________________________ > >> > grass-dev mailing list > >> > [email protected] > >> > https://lists.osgeo.org/mailman/listinfo/grass-dev --000000000000ea919305b9a91aab Content-Type: text/html; charset="UTF-8" Content-Transfer-Encoding: quoted-printable <div dir=3D"ltr"><div>A note about the usage of GRASS_EPSILON =3D 1.0e-15:<= br><br>The range of single precision floating point is about 1.4012984643e= =E2=88=9245 to 3.4028234664e38<br></div><div>That means for large numbers, = using GRASS_EPSILON will also detect floating point precision limits, not o= nly meaningful differences. For small numbers, GRASS_EPSILON will not detec= t meaningful differences. GRASS_EPSILON could instead be modified with some= thing like "max(abs(map_A), abs(map_B)) * 1.0e-15" to test for me= aningful differences.</div><div><br></div><div>Markus M<br></div><div><br><= /div><div><br></div><div>On Sun, Jan 24, 2021 at 5:32 PM Markus Metz <<a= href=3D"mailto:[email protected]">[email protected]= m</a>> wrote:<br>><br>><br>><br>> On Sun, Jan 24, 2021 at 5:= 06 PM ming han <<a href=3D"mailto:[email protected]">[email protected]= </a>> wrote:<br>> ><br>> > Hi Everyone<br>> ><br>> = > =C2=A0 =C2=A0Many thanks for your help. Is if(fabs(map_A - map_B) <= =3D 1.0e-15, ... ) approach may increase the runtime compare to '=3D=3D= ' way?<br>><br>> The formulas are very simple, I don't think = that differences in runtime can be reliably measured. (De-)compression of t= he data and the operating system's file cache have a much stronger infl= uence on the runtime.<br>><br>> Markus M<br>> ><br>> > Th= anks<br>> > Ming<br>> ><br>> > Markus Metz <<a href=3D= "mailto:[email protected]">[email protected]</a>>= ; =E4=BA=8E2021=E5=B9=B41=E6=9C=8824=E6=97=A5=E5=91=A8=E6=97=A5 =E4=B8=8A= =E5=8D=8810:57=E5=86=99=E9=81=93=EF=BC=9A<br>> >><br>> >>= Trying to answer the original question: with a DCELL map "cat1_acc_ri= v" and a FCELL map "cat1_minacc", why is "float(cat1_ac= c_riv) =3D=3D float(cat1_minacc)" not equal to "int(cat1_acc_riv)= =3D=3D int(cat1_minacc)" ?<br>> >><br>> >> int trun= cates to integer while float converts to single precision floating point. E= .g. with cat1_acc_riv =3D 1.1 and cat1_minacc =3D 1.9, "float(cat1_acc= _riv) =3D=3D float(cat1_minacc)" becomes "1.1 =3D=3D 1.9" wh= ereas "int(cat1_acc_riv) =3D=3D int(cat1_minacc)" becomes "1= =3D=3D 1", thus the results are different.<br>> >><br>> &= gt;> Another reason for possible differences is that float can only repr= esent max 7 decimal digits. E.g. float(194320567) becomes 194320560 but int= (194320567) preserves the value 194320567.<br>> >><br>> >>= ; Thus the safest is to cast everything to the type with the highest precis= ion. In this case with FCELL and DCELL, use "double(cat1_acc_riv) =3D= =3D double(cat1_minacc)" or even better the suggestion of Markus N.<br= >> >><br>> >> Markus M<br>> >><br>> >><= br>> >> On Sun, Jan 24, 2021 at 3:51 PM ming han <<a href=3D"ma= ilto:[email protected]">[email protected]</a>> wrote:<br>> >>= ><br>> >> > Hi Markus and Micha<br>> >> ><br>&g= t; >> > =C2=A0 =C2=A0 =C2=A0I am just trying to find grids have th= e same values in these two rasters, I will try the threshold approach.<br>&= gt; >> ><br>> >> > Thanks<br>> >> > Ming<b= r>> >> ><br>> >> > Markus Neteler <<a href=3D"ma= ilto:[email protected]">[email protected]</a>> =E4=BA=8E2021=E5=B9=B41= =E6=9C=8824=E6=97=A5=E5=91=A8=E6=97=A5 =E4=B8=8A=E5=8D=886:58=E5=86=99=E9= =81=93=EF=BC=9A<br>> >> >><br>> >> >> Hi Ming= ,<br>> >> >><br>> >> >> On Sun, Jan 24, 2021 = at 10:49 AM ming han <<a href=3D"mailto:[email protected]">dustming@gma= il.com</a>> wrote:<br>> >> >> ><br>> >> >&= gt; > Hi Micha<br>> >> >> ><br>> >> >> = > =C2=A0 =C2=A0 =C2=A0Many thanks for your reply.<br>> >> >&= gt; > =C2=A0 =C2=A0 =C2=A0Here is the command I am using:<br>> >&g= t; >> ><br>> >> >> > =C2=A0 =C2=A0 =C2=A0if(floa= t(cat1_acc_riv) =3D=3D float(cat1_minacc), str_r, null())<br>> >> = >> ><br>> >> >> > =C2=A0 =C2=A0 =C2=A0 The str_r= is a CELL raster. the result is different when I change it to:<br>> >= ;> >> > =C2=A0 =C2=A0 =C2=A0 =C2=A0if(int(cat1_acc_riv) =3D=3D = int(cat1_minacc), str_r, null())<br>> >> >><br>> >>= >> Note that numerical "equality" is better tested with a = threshold test<br>> >> >> against the map pixel difference.<= br>> >> >> As the threshold, we use GRASS_EPSILON which is d= efined as 1.0e-15.<br>> >> >><br>> >> >> Henc= e the test needs to be implemented in a different way, i.e. by<br>> >= > >> using an epsilon.<br>> >> >> Essentially somet= hing like this:<br>> >> >><br>> >> >> if(fabs= (map_A - map_B) <=3D 1.0e-15, ... )<br>> >> >><br>> &g= t;> >> In your case (untested):<br>> >> >> r.mapcal= c diffepsilon =3D if( abs( map_A - map_B) <=3D 1.0e-15, str_r , null())<= br>> >> >><br>> >> >> See related discussions= here: [1], [2] and elsewhere.<br>> >> >><br>> >> &= gt;> [1] Comment by Glynn: <a href=3D"https://trac.osgeo.org/grass/ticke= t/2854#comment:9">https://trac.osgeo.org/grass/ticket/2854#comment:9</a><br= >> >> >> [2] Comment by Glynn:<br>> >> >> <a = href=3D"https://lists.osgeo.org/pipermail/grass-user/2015-October/073200.ht= ml">https://lists.osgeo.org/pipermail/grass-user/2015-October/073200.html</= a><br>> >> >><br>> >> >> Best,<br>> >&g= t; >> Markus<br>> >> ><br>> >> > ____________= ___________________________________<br>> >> > grass-dev mailing= list<br>> >> > <a href=3D"mailto:[email protected]">gr= [email protected]</a><br>> >> > <a href=3D"https://lists.= osgeo.org/mailman/listinfo/grass-dev">https://lists.osgeo.org/mailman/listi= nfo/grass-dev</a></div></div> --000000000000ea919305b9a91aab-- --===============6778851474218773543== Content-Type: text/plain; charset="us-ascii" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit Content-Disposition: inline _______________________________________________ grass-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/grass-dev --===============6778851474218773543==--