[alg] snelle asin,acos,atan.

[Alarmnummer]

Ik ben op dit moment bezig met j2me (java 2 microedition) bezig om een spel te schrijven voor pda`s en snelle mobiele telefoons. j2me heeft geen beschikking over floats en je zult zelf een oplossing hiervoor vinden.

Op dit moment maak ik gebruik van een fixedpoint library waar ik een float kan coderen in een integer. Op zich werkt het heel aardig, maar ik zit 30% van de tijd in de atan functie en dit vind ik te veel van het goeie.

Een sin,cos en tan kan ik wel snel maken door te werken met een lookup table. Floating point getallen wil ik verder vermijden door te werken in kleinere eenheden (10e van graden ipv graden, cm ipv m). Maar een goeie oplossing voor de atan heb ik nog niet kunnen vinden.

Ik heb net .oisyn en soultaker gevraagd en ze kwamen ongeveer op hetzelfde uit. Dmv een search in een tree op zoek te gaan naar de hoek die het beste past bij een bepaalde x en y. Helaas ben je dan wel een aantal stappen nodig om de bepalen welke hoek je hebt. Mijn vraag is of mensen nog een betere manier weten.

[Macros]

atan is toch een geinverteerd tan? Kan je niet een lookup table maken in een switch statement? Zo had iemand al zijn lookuptables geconverteerd naar switches, omdat hij 10 x zoveel code memory had dan data memory.

[.oisyn]

Waarom kun je de lookuptable methode eigenlijk niet gebruiken voor acos, asin en atan?

[Alarmnummer]

.oisyn schreef op 19 juli 2004 @ 17:28:
Waarom kun je de lookuptable methode eigenlijk niet gebruiken voor acos, asin en atan?

Als je werkt met bv 360 graden (of 720 halve graden) dan kan je eenvoudig het volgende doen:

code:

int tan(int hoek){
      return _tanTable[hoek%360];
}

Maar hoe wou je de verhouding tussen een x en een y gaan normaliseren?

code:

int atan(int x, int y){
      return _atanTable[?];
}

[ Voor 9% gewijzigd door Alarmnummer op 19-07-2004 17:31 ]

[.oisyn]

Da's atan2, en die maakt weer gebruik van atan door y/x te nemen, en bovendien te controleren of x == 0 en in welk kwadrant het resultaat terecht zou moeten komen door x en y te testen tegen 0

.edit: ik zou je hoeken trouwens definieren in modulo een macht van 2, zodat je de & kunt gebruiken ipv de vrij dure %

[ Voor 25% gewijzigd door .oisyn op 19-07-2004 17:34 ]

[Verwijderd]

Ik weet niet of ik op het goeie spoor zit (of ik de vraag juist heb gelezen...), maar toch even proberen
atan (tegenwoordig spreken we wel over bgtan heeft als domein gans R, maar het beeld ligt tss -pi/2 en +pi/2. Er zijn 2 horizontale asymptoten, nl de gegeven randelementen van het beeld.
De grootste functieschommelingen liggen in het bereik -900° , + 900° (ff mijn interpretatie).
Dus in dit bereik zou je eventueel een tabelletje kunnen aanleggen van hoeken.
Buiten dit bereik zou ik opteren om het exponentieel aan te pakken, 1000°, 10000°, 1000000,...

En ik zie niet in waarom je dit zegt: int atan(int x, int y){
Je kan toch makkelijk de hoek berekenen waarom het gaat?

[Macros]

Ik denk niet meer in graden, alles wat ik met deze zooi krijg is altijd gewoon in radialen (zoals het hoort).
Atan is een simpele functie. Alleen in het midden is het interesant. Als je aan de randen komt, dus bv. < -20 en > 20, dan is het gewoon -1 respectievelijk +1.
Als je een switch gebruikt kun je doen (ik gebruik hier even een switch op een float, kan je zelf wel omzetten naar fixed point int):

Java:

switch(x){
   default: // groter dan 20 ofzo
      return 1;
   case 0:  // midden van de functie
      return 0;
   case 0.1: // stijgende lijn, ongeveer helling 1
      return 0.1;
   case 0.2:
      return 0.19;
etc....
}

In je switch, switch je alleen op positieve getallen. De functie die je aanroept, moet zelf maar zorgen dat het antwoord omzet als het toevallig een negatieve x is (scheelt de helft qua grote van switch).

[Verwijderd]

Macros schreef op 19 juli 2004 @ 17:45:
Ik denk niet meer in graden, alles wat ik met deze zooi krijg is altijd gewoon in radialen (zoals het hoort).
Atan is een simpele functie. Alleen in het midden is het interesant. Als je aan de randen komt, dus bv. < -20 en > 20, dan is het gewoon -1 respectievelijk +1.
Als je een switch gebruikt kun je doen (ik gebruik hier even een switch op een float, kan je zelf wel omzetten naar fixed point int):

Java:

switch(x){
   default: // groter dan 20 ofzo
      return 1;
   case 0:  // midden van de functie
      return 0;
   case 0.1: // stijgende lijn, ongeveer helling 1
      return 0.1;
   case 0.2:
      return 0.19;
etc....
}

Jij denkt aan er op dezelfde manier over als ik (ik had het trouwens voor de makkelijkheid ff trg omgezet naar graden).

Alarmnummer: in je andere topics staat ergens dat je werkt met een soort van orthonormaal assenstelsel, dus zoals Macros zegt, kunnen we dan de helft van de berkening/tabel scrhappen door gewoon ff te checken als volgt:

hoek==berekeningetje om positieve hoek te bekomen mbv van x en y
If(y>0){
bgtan==zoekbgtanintabel(hoek);
}
else
{
bgtan==-zoekbgtanintabel(hoek);
};

mm dat was dus java in mn uber-1337 javaskills (moet dringend eens java leren)

Door ff een snel if then else-etje te doen, heb je nog maar de helft van de tabelgrootte nodig.

[ Voor 1% gewijzigd door Verwijderd op 19-07-2004 18:06 . Reden: /quote vergeten ]

[Soultaker]

Zoals gezegd is tan een repeterende functie met een periode pi. Het is bovendien een 'oneven' functie, wat betekent dat tan(x) = -tan(-x). Alleen het functiebereik van 0 tot pi/2 is dus relevant.

Daar is niet zo makkelijk een switch tabel van te maken doordat de distributie tussen invoer- en uitvoerwaarden nogal verschillend is. Je kunt echter wel een tabel maken die ofwel de invoer, ofwel de uitvoer een vaste precisie geeft. Voor de functie tan(x) kun je bijvoorbeeld de invoerwaarden steeds een vast verschil geven.

x	0	pi/8	pi/4	3*pi/8	pi/2
tan(x)	0	0,3926	0,4142	2,4142	+INF

Om atan te berekenen kan je die tabel omgekeerd gebruiken: als het de invoerwaarde y dicht bij 0,41 in de buurt zit, dan zal atan(y) dicht bij pi/4 in de buurt zitten. Omdat de waarden in de tweede kolom, waar je in zoekt, oplopend gesorteerd zijn kun je er met een binary search in zoeken.

In Java-code ziet dat er bijvoorbeeld zo uit:

Java:

public class Atan
{
    static final int samples = 256;
    static final float limits[] = new float[samples];
    static final float values[] = new float[samples];
    
    static {
        for(int n = 0; n < samples; ++n)
        {
            limits[n] = (float)Math.tan( ((n + 1.0) * (Math.PI/2.0)) / (samples + 1) );
            values[n] = (float)((n + 0.5) * (Math.PI/2.0)) / (samples + 1);
        }
    }
    
    public static float atan(float x)
    {
        if(x < 0)
            return -atan(-x);
        
        if(x < limits[0])
            return 0;

        if(x >= limits[samples-1])
            return (float)(Math.PI/2.0);

        int lo = 0, hi = samples - 1;
        while (lo < hi)
        {
            int mid = (lo + hi)/2;
            if(x >= limits[mid])
                lo = mid + 1;
            else
                hi = mid;
        }
        return values[lo];
    }
}

De nauwkeurigheid van de functie is ten dele afhankelijk van de manier waarop de return value berekend wordt en ten dele (natuurlijk) van het aantal samples. Bij deze code is de afwijking tussen de uitvoer van deze implementatie en de werkelijke waarde van atan(x) maximaal (pi/2)/samples; voor 256 samples is dit dus ongeveer 0,0061. De gemiddelde afwijking kan verder omlaag worden gebracht door (bijvoorbeeld) het resultaat lineair te interpoleren tussen twee samples. Het is de vraag of dat de extra code waard is, natuurlijk.

[Soultaker]

Macros schreef op 19 juli 2004 @ 17:45:
Atan is een simpele functie. Alleen in het midden is het interesant. Als je aan de randen komt, dus bv. < -20 en > 20, dan is het gewoon -1 respectievelijk +1.
Als je een switch gebruikt kun je doen (ik gebruik hier even een switch op een float, kan je zelf wel omzetten naar fixed point int):

Java:

switch(x){
   default: // groter dan 20 ofzo
      return 1;
   case 0:  // midden van de functie
      return 0;
   case 0.1: // stijgende lijn, ongeveer helling 1
      return 0.1;
   case 0.2:
      return 0.19;
etc....
}

Dit is (in het gunstigste geval; dwz. als Java er een jumptable voor gebruikt) toch ook gewoon een constante lookup table?

[Verwijderd]

Soultaker schreef op 19 juli 2004 @ 18:11:
Zoals gezegd is atan een repeterende functie met een periode pi. Het is bovendien een 'oneven' functie, wat betekent dat atan(x) = -atan(-x). Alleen het functiebereik van 0 tot pi/2 is dus relevant.

de Tangensfunctie is een repeterende functie, de boogtangensfunctie niet, het is de inversie functie van de beperkte tangensfunctie tot het bereik -pi/2, pi/2.

(*aspfreakout schopt soultaker een paar jaartjes trg naar school )
meer info: http://www.wisfaq.nl/fram.../showrecord3.asp?id=23331

[Soultaker]

Ik had het in die alinea dus over de tangens; die a'tjes hoorde daar gewoon niet. (De tabel die ik gaf was ook van de tanges en die klopt wel.)

[Soultaker]

Oh trouwens, ik bedenk me net dat atan(x) = pi/2 - atan(1/x) als x>=1. Op deze manier kun je je invoer dus eenvoudig beperken tot een bereik van 0 tot 1 en daar is weer uitstekend een constante lookup table voor te maken, denk ik?

[Verwijderd]

idd, maar volgens mij gaat dit ongeveer op hetzelfde neerkomen want je gaat een grotere nauwkeurigheid nodig hebben door 1/x.
Maar dit lijkt wel de beste aanpak omdat je anders nog ff een bereik moet gaan bedenken waarin je een table gaat gebruiken. Nu gaan alle waarden ongeveer even onnauwkeurig zijn, met een zelfgekozen bereik is de nauwkeurigheid niet direct voorspelbaar...

[WVL_KsZeN]

i quote even een stukje code van een deel van een spelletje van me. Dit is geschreven in 6502 assembler, maar ik hoop dat het je duidelijk is wat er gebeurt

code:

         ;calculate angle in which ball is travelling
         ;using angle=atan(y/x) -> y/x is done using logarithms log(y/x)=log(y)-log(x)
         ;the atan-table has a inverse log built in, to undo log(y/x)

         ldy spd2ylow
         ldx spd2xlow

         ;if x=0 then angle = $40, as log(0) doesn't exist (== -infinity)

         bne *+15
         cpy #0          ;if BOTH speeds are zero
         bne *+7         ;don't calculate the friction
         tya
         sta angle
         sec             ;skip adding friction
         rts
         lda #$40
         bne undoswap

         lda loglow,y   ;calculate log(y/x)
         sec
         sbc loglow,x
         lda loghigh,y
         sbc loghigh,x  ;use only resulting loghigh (y/x) for calculating angle
         tay
         lda atandata,y ;do angle=atan(invlog(log(y/x)

dit is zo'n beetje de simpelste nauwkeurige methode met een kleine look-up table...

suc6!

edit:
wat meer uitleg : ik comprimeer de lookup tabel van de atan dus door een log en inverse log toe te passen, op die manier krijg ik ook nog eens gratis een x/y cadeau. Ik kan me niet meer precies herinneren welke log en inverse log het precies is, iig niet log2 of log10, maar ik heb een beetje gepuzzeld zodat de hele tabel precies van 0 tot 255 liep...

[ Voor 15% gewijzigd door WVL_KsZeN op 19-07-2004 19:49 ]

[Verwijderd]

Is het niet mogelijk om de taylor functie van atan te gebruiken ??

atan(x) = x - x^3/3 + x^5/5 - x^7/7 + x^9/9 ...........

[Cavorka]

Emptor22: Werkt die reeks van jou niet alleen om nul?

Chebychev polynomen misschien? Zo werken rekenmachines in ieder geval wel...

[Sjeik]

Verwijderd schreef op 19 juli 2004 @ 19:49:
Is het niet mogelijk om de taylor functie van atan te gebruiken ??

atan(x) = x - x^3/3 + x^5/5 - x^7/7 + x^9/9 ...........

Dit is idd de oplossing lijkt mij. Deze taylorreeksen gaan supersnel naar een heel nauwkeurige benadering van het getal dat je wil hebben.

[Macros]

Taylor is ook maar een benadering en niet erg efficient. Hij moet gewoon een lookup table maken, en dan zelf bepalen hoe nauwkeurig hij die wilt hebben

[Verwijderd]

Arctg(x) = sigma n [2(wortel2-1)2n+1 (-1)n / (2n+1)] T2n+1(x)
= 2(Ö2-1) sigma n [(2wortel2-3)n / (2n+1)] T2n+1(x)
(heb het mr even opgezocht
Lijkt mij toch niet handig om dit telkens te laten berekenen.

Nu ja, Allarmnummer zegt dat hij het wil tot op een tiende van een graad. Aangezien 0.1 ongeveer 0.0017 raden is, heeft hij een tabel met 570 waarden nodig om alle bgtangensen te berekenen die hij nodig heeft!

[MSalters]

Hoe dat zo? Volgens mij heb je het bereik [0, pi/4] nodig (zoals Soultaker terecht opmerkt) en dat is 45 graden of 450 waarden in een lookup. Binary search gebruiken, Newton is ook mogelijk maar biedt relatief weinig voordelen hier. Het verschil is ~9 stappen om ~6 stappen, en zonder floats is Newton per stap makkelijk 50% duurder.

[SuperRembo]

Het probleem is dat je niet atan(x) zoekt, maar atan(y/x). En dat gaat niet 1 2 3 in een lookup.

[Soultaker]

Waarom is dat fundamenteel lastiger dan? Ok, je moet controleren of x misschien 0 is en de deling daadwerkelijk uitvoeren, maar daarna kun je alsnog een lookup doen. (Overigens moet je nog net iets meer doen om het gedrag van atan2() te krijgen, maar dat is allemaal verwaarloosbaar).

[SuperRembo]

Hier wordt een fixed point benadering gegeven (voor x,y>0)

code:

r = (x-y)/(x+y)
//benadering 1, afwijking max 0.07 rad 
theta1 = (1 - r) * pi/4
//benadering 2, afwijking max 0.01 rad
theta2 = 0.1963 * r^3 - 0.9817 * r + pi/4

De afleiding is me niet echt duidelijk.

[.oisyn]

Een oplossing voor atan2 zou ook kunnen door je 2d ruimte steeds door 2en te delen adhv een lijn zodat je uiteindelijk (oneindig grote) taartstukken opruimt. Elk taartstuk representeert een bepaald antwoord, het opzoeken is dan simpelweg je x en y coördinaat door de ontstane bsp boom te halen. Deze O(log n) oplossing is natuurlijk wat duurder dan de O(1) voor een simpele lookup, maar erg veel checks zijn er natuurlijk niet nodig (360 individuele graden kost je slechts 9 checks)

[Zoijar]

Deze heb ik ook nog gevonden:

code:

r(x) = (x + 0.0.43157974*x^3)/(1 + 0.76443945*x^2 + 0.05831938*x^4)

For -1 <= x <= 1 we have
|arctan(x) - r(x)| < 0.00000642

.oisyn schreef op 20 juli 2004 @ 00:24:
Een oplossing voor atan2 zou ook kunnen door je 2d ruimte steeds door 2en te delen adhv een lijn zodat je uiteindelijk (oneindig grote) taartstukken opruimt. Elk taartstuk representeert een bepaald antwoord, het opzoeken is dan simpelweg je x en y coördinaat door de ontstane bsp boom te halen. Deze O(log n) oplossing is natuurlijk wat duurder dan de O(1) voor een simpele lookup, maar erg veel checks zijn er natuurlijk niet nodig (360 individuele graden kost je slechts 9 checks)

Onvoorspelbare conditionals, en memory accesses kunnen dat heel erg traag maken. Een branch prediction failure (die je zeker vaak krijgt) of cache misses kosten aanzienlijk tijd.
(hoewel ik niet precies weet hoe dat zit voor mobiele telefoons. Maar veel pda's draaien gewoon op intel pentiums toch?)

[ Voor 76% gewijzigd door Zoijar op 20-07-2004 11:04 ]

[Alarmnummer]

Zoijar schreef op 20 juli 2004 @ 10:53:
(hoewel ik niet precies weet hoe dat zit voor mobiele telefoons. Maar veel pda's draaien gewoon op intel pentiums toch?)

PDA gebruiken geen pentiums, meer XScale/Arm processoren.

[.oisyn]

Zoijar schreef op 20 juli 2004 @ 10:53:
Onvoorspelbare conditionals, en memory accesses kunnen dat heel erg traag maken. Een branch prediction failure (die je zeker vaak krijgt) of cache misses kosten aanzienlijk tijd.
(hoewel ik niet precies weet hoe dat zit voor mobiele telefoons. Maar veel pda's draaien gewoon op intel pentiums toch?)

Je hebt het hier over java code op mobile devices, dan heb je daar niet echt veel last van omdat het gewoon heel erg schaalbare processors zijn, en het jitten sowieso al veel overhead met zich meebrengt. Bovendien kan het ook zonder conditionals door de checks achter elkaar te doen en de antwoorden te coderen in de bits van een int. Het uiteindelijke antwoord kun je dan gewoon opzoeken in een array.