2011-08-16

Supersampling

Supersamling polega na tym by na każdy piksel obrazu wygenerować więcej subpikseli niż nam potrzeba, a następnie uśrednić ich wartość. Uśrednianie to proces resamplingu i istnieje wiele metod, a w zasadzie filtrów resamplujących. W najprostszej postaci bierzemy średnią ze wszystkich subpikseli.

Supersampling + resampling pozwala nam na polepszenie jakości sceny. Usuwamy z niej schodki na figurach, schodki na teksturach, interferencje na teksturach, które manifestują się jeśli tekstura jest jakimś regularnym wzorem.

Przypuśćmy, że mamy sygnał dźwiękowy próbkowany z częstotliwością 44KHz (przez analogię możemy go traktować jako linię obrazu). Chcemy go przekształcić w sygnał próbkowany z częstotliwością 14KHz (przez analogię nasz element tekstury zajmuje niewielki obszar ekranu, trzeba go mniej więcej 3x pomniejszyć). Uwaga próbkowanie z częstotliwością 44KHz oznacza, że sygnał niesie informację o częstotliwości 22KHz i mniejszych. Jeśli wybierzemy co mniej więcej 3 próbkę dźwięk będzie zniekształcony. Do sygnału wyjściowego wkradną się składowe widma o częstotliwościach powyżej 7KHz. Dokładniej zostaną one odbite tak, że 8KHz pojawi się na 6KHz, 16KHz pojawi się na 5KHz. I już nic z tym nie zrobimy, nie jesteśmy w stanie odróżnić prawdziwego sygnału 5KHz od tego odbitego. Aby tego uniknąć przed downsamplingiem sygnał filtruje się filtrem dolnoprzepustowym o częstotliwości 7KHz. Tak robi się w przypadku dźwięku i filtrów resamplujących typu Lanczos i Bicubic, które są pewnym przybliżeniem filtru dolnoprzepustowego. Generalnie w przypadku dźwięku wszelkie przybliżenia są niepożądane, gdyż ucho ludzkie w przeciwieństwie do oka analizuje widmo sygnału. W przypadku obrazu możemy wziąć 3 próbki i je uśrednić (to ciągle jest swego rodzaju dosyć niedoskonały filtr dolnoprzepustowy), w przypadku dźwięku nie.

Czyli myśląc o wygładzaniu, filtrowaniu, samplowaniu, resamplowaniu powinniśmy raczej myśleć o częstotliwościach. To, że zwykłe uśrednianie też działa (ale nie daje dobrego rezultatu) wynika z faktu jak działa oko.

Mamy sinusoidę powiedzmy o f=1KHz, próbkowanie fs=44KHz. Nasza sinusoida ewidentnie ma schodki. Ale do reprodukcji tych schodków potrzeba częstotliwości wyższych niż 44KHz, a tych nie w naszym sygnale. Jeśli taki schodek z DAC puścimy na głośniku to będziemy mieć zniekształcenia, ale jeśli użyjemy filtru dolnoprzepustowego 22KHZ to na wyjściu zostanie w idealnym świecie pojedyncza sinusoida (pik w domenie częstotliwości).

Innymi słowy filtr dolnoprzepustowy wygładził nam schodki usuwając z sygnału wyższe częstotliwości. Dokładnie to samo robimy w resamplingu.

By móc bez zniekształceń przedstawić sygnał o wyższej częstotliwości w postaci sygnału o niższej częstotliwości musimy mieć więcej próbek niż potrzebujemy.

Idealnie chcąc uniknąć jakichkolwiek błędów potrzebujemy po próbkę na każdy piksel tekstury inaczej zawsze już podczas supersamplingu wprowadzamy do sygnału zakłócenia (aliasing). Niestety z uwagi na złożoność obliczeniową jest to dosyć kłopotliwe. Dlatego np. przyjmując że samplujemy 9x9 subpikseli na piksel nasz obraz wzorcowy posiada już w sobie zakłócenia. Proces resamplingu nie usuwa ich. Pozwala nam się pozbyć tych związanych z częstotliwościami które są w 9x9, a nie powinno ich być w 1x1. W zależności od odległości tekstury od ekranu usuwamy tym tylko pewną część zakłóceń, ale można powiedzieć tą najważniejszą znacząco poprawiającą obraz.

Supersampling działa kompleksowo na całą scenę. Na tekstury (wygładzanie + interferencje), na brzegi obiektów (wygładzenie), na wiele małych obiektów równo rozłożonych.

Z uwagi na to, że nie wszystkie piksele w obrazie końcowym wymagają takiego samego supersamplingu stosuje się samplowanie adaptacyjne, w którym uwzględnia się różnice pomiędzy sąsiednimi pikselami i jeśli jest ona duża to obszar dzieli się na mniejsze subpiksele. To problem na zupełnie inny wpis i czas (implementacji).

By walczyć z zakłóceniami wprowadzanymi podczas tworzenia obrazu subpikseli stosuje się samplowanie losowe. Staramy się nasze promienie wysyłać nieregularnie rozłożone.

Jeśli promienie są generowane regularnie w każdy subpiksel, to mamy tzw. samplowanie typu Grid, jeśli wysyłamy je losowo, ale ciągle w granicach subpiksela to mamy samplowanie typu Jitter.

Grid:


  
  
  
    
      
        image/svg+xml
        
        
      
    
  
  
    
    
    
    
    
    
    
    
    
    
    
    
    
  

Jitter:


  
  
  
    
      
        image/svg+xml
        
        
      
    
  
  
    
    
    
    
    
    
    
    
    
    
    
    
    
  

Zupełnie inną klasą samplerów są te które nie wysyłają po jednym promieniu na subpiksel. Mówimy wtedy o samplowanie niejednorodnym typu Random (zupełnie losowe), Poisson Disk (losowe, ale równomiernie rozmieszczone, podobnie jak czopki w oku).

Uwaga supersampling+resampling, a filtrowanie tekstur (też wygładzanie) to zupełnie co innego. W filtrowaniu tekstur wyliczamy wartości koloru w konkretnej współrzędnej zmiennoprzecinkowej wykorzystując wartości w punkach całkowitych (lepiej o punktach myśleć jak o kwadracikach o takim samym kolorze). Informacje którą generujemy tam nie ma. Zgadujemy co znajduje się w punktach pośrednich. Inaczej niż do supersamplingu gdzie wartości subpikseli są wyliczane. Filtrowanie wygładzi nam teksturę ale nie usunie interferencji i schodków na brzegach obiektów.

Poniżej kod samplerów typu Jitter i Grid. Innych na razie nie implementowałem. W przyszłości może pojawi się samplowanie adaptywne. Niejednorodne to śpiew przyszłości.

Kod:
public enum SampleMethod
{
    Grid,
    Jitter
}

internal class PixelInfo
{
    /// <summary>
    /// Within RenderOptions.Width and RenderOptions.Height.
    /// </summary>
    public readonly Vector2 RayPos;

    /// <summary>
    /// Within RenderOptions.Width and RenderOptions.Height times sampling size.
    /// </summary>
    public readonly Vector2 PixelPos;

    public PixelInfo(Vector2 a_ray_pos, Vector2 a_pixel_pos)
    {
        Debug.Assert(a_pixel_pos.X.Fraction().IsZero());
        Debug.Assert(a_pixel_pos.Y.Fraction().IsZero());

        RayPos = a_ray_pos;
        PixelPos = a_pixel_pos;
    }
}

internal abstract class Sampler
{
    private RenderOptions m_render_options;

    public static Sampler Create(RenderOptions a_render_options)
    {
        if (a_render_options.SampleMethod == SampleMethod.Grid)
            return new GridSampler(a_render_options);
        if (a_render_options.SampleMethod == SampleMethod.Jitter)
            return new JitterSampler(a_render_options);
        else
            throw new InvalidOperationException();
    }

    public Sampler(RenderOptions a_render_options)
    {
        m_render_options = a_render_options;
    }

    public IEnumerable<PixelInfo> GetSamples()
    {
        return GetSamples(new Rectangle(0, 0, m_render_options.Width,
            m_render_options.Height));
    }
    public IEnumerable<PixelInfo> GetSamples(Rectangle a_rect)
    {
        for (int y = a_rect.Top; y < a_rect.Bottom; y++)
        {
            for (int yy = 0; yy < m_render_options.SampleSize; yy++)
            {
                for (int x = a_rect.Left; x < a_rect.Right; x++)
                {
                    for (int xx = 0; xx < m_render_options.SampleSize; xx++)
                    {
                        yield return CreatePixelInfo(x, y, xx, yy,
                            m_render_options.SampleSize);
                    }
                }
            }
        }
    }

    protected abstract PixelInfo CreatePixelInfo(int a_x, int a_y, 
        int a_xx, int a_yy, int a_size);
}

internal class JitterSampler : Sampler
{
    private Random m_random;

    public JitterSampler(RenderOptions a_render_options)
        : base(a_render_options)
    {
        m_random = new Random(a_render_options.RandomSeed);
    }

    protected override PixelInfo CreatePixelInfo(int a_x, int a_y, 
        int a_xx, int a_yy, int a_size)
    {
        return new PixelInfo(
            new Vector2(
                a_x + (a_xx + m_random.NextDouble()) / a_size,
                a_y + (a_yy + m_random.NextDouble()) / a_size),
            new Vector2(a_x * a_size + a_xx, a_y * a_size + a_yy)
        );
    }
}

internal class GridSampler : Sampler
{
    public GridSampler(RenderOptions a_render_options)
        : base(a_render_options)
    {
    }

    protected override PixelInfo CreatePixelInfo(int a_x, int a_y, 
        int a_xx, int a_yy, int a_size)
    {
        return new PixelInfo(
            new Vector2(
                a_x + (a_xx + 0.5) / a_size,
                a_y + (a_yy + 0.5) / a_size),
            new Vector2(a_x * a_size + a_xx, a_y * a_size + a_yy)
        );
    }
}
Przykłady samplowania z różną dokładnością. Resamplowanie zawsze typu Box. Zauważmy, że Jitter nadaję się tylko dla większej ilości próbek, gdyż inaczej brzegi są poszarpane. Widzimy, że czy większa ilość próbek tym zakłócenia pojawiają się coraz bliżej horyzontu, w momencie w którym częstotliwość próbkowania zrównuje się z częstotliwością biało-czerwonej siatki.


Grid 1x1:



Jitter 1x1:



Grid 3x3:



Jitter 3x3:



Grid 6x6:



Jitter 6x6:



Grid 9x9:



Jitter 9x9:



2011-08-12

Pobranie wszystkich metod typu

Wszystkich czyli publicznych, chronionych, prywatnych. Także tych ukrytych. Pomijamy metody abstrakcyjne. Z nadpisanym metod wirtualnych zwracamy tylko te "najmłodsze" w hierarchii klas.

Możemy także zdecydować się czy odfiltrować metody generowane dla automatycznych właściwości.

Metody są zwracane w sposób posortowany najpierw po typie (od najmłodszego do bazowego), a dalej w kolejności występowania.

I znowu problemy. Tym razem przypuśćmy, że mamy metodę publiczną zdefiniowaną w klasie A, mamy klasę B wywodzącą się z A. Pobieramy tą metodę z typu A i z typu B, chcemy ustalić, że są takie same. Ich porównanie nic nie da. Ich ReflectedType są różne. Pozostaje więc ręcznie wykazać, że metody mają takie same sygnatury. Z całego dobrodziejswta MethodInfo decydujemy się nazwy i parametry. Niestety, okaże się parametrów też nie możemy porównać. Żeby było śmieszniej różnią się tylko MetadataToken. Z całego dobrodziejstwa ParameterInfo decydujemy się na porównanie ParameterType. Bez problemu radzimy sobie z rozróżnieniem parametrów z out, in, ref. Zgodnie z moimi testami wszystko działa. Ale nie wykluczam, że z uwagi na to, że do porównania nie biorę wszystkich elementów istnieją takie których porównanie zwróci fałszywy wynik.

Kod:
public static IEnumerable<type> GetBaseTypes(this Type a_type, 
    bool a_with_this = false)
{
    if (a_with_this)
        yield return a_type;

    Type t = a_type;

    while (t.BaseType != null)
    {
        t = t.BaseType;
        yield return t;
    }
}

public static IEnumerable<methodinfo> GetBaseDefinitions(this MethodInfo a_mi, 
    bool a_with_this = false)
{
    if (a_with_this)
        yield return a_mi;

    MethodInfo t = a_mi;

    while ((t.GetBaseDefinition() != null) && (t.GetBaseDefinition() != t))
    {
        t = t.GetBaseDefinition();
        yield return t;
    }
}

public static bool IsDerivedFrom(this MethodInfo a_mi, MethodInfo a_base,
    bool a_with_this = false)
{
    if (a_mi.Name != a_base.Name)
        return false;
    if (a_mi.DeclaringType == a_base.DeclaringType)
    {
        if (!a_mi.GetParameters().Select(p => p.ParameterType).SequenceEqual(
            a_base.GetParameters().Select(p => p.ParameterType)))
        {
            return false;
        }

        return a_with_this;
    }

    return a_mi.GetBaseDefinitions().Contains(a_base);
}

public static IEnumerable<MethodInfo> GetAllMethods(this Type a_type,
    bool a_include_autogenerated = false)
{
    List<MethodInfo> result = new List<MethodInfo>();

    foreach (var t in a_type.GetBaseTypes(true))
    {
        if (t == typeof(Object))
            break;
        if (t == typeof(ValueType))
            break;

        MethodInfo[] methods = t.GetMethods(
            BindingFlags.Public | BindingFlags.NonPublic |
            BindingFlags.Instance);

        foreach (var method in methods.Reverse())
        {
            if (method.IsDefined(typeof(CompilerGeneratedAttribute), true) &&
                !a_include_autogenerated)
            {
                continue;
            }

            if (method.IsAbstract)
                continue;
            if (result.All(m => !m.IsDerivedFrom(method, true)))
            {
                result.Add(method);
                yield return method;
            }
        }
    }
}

2011-08-11

Pobranie wszystkich właściwości i pól typu

Wszystkich czyli także tych ukrytych, prywatnych, z typu i wszystkich podtypów i co ważniejsze unikalnych. W przypadku właściwości nadpisanych wirtualnych zwracamy tylko jedną. Pomijamy właściwości abstrakcyjne.

Zarówno właściwości jak i pola staramy się zwrócić posortowane według hierarchii typów, a następnie posortowane według tego jak występują w Typie, co odzwierciedla ich porządek w kodzie.

W przypadku pól istnieje możliwość odfiltrowania pól dodanych przez właściwości automatyczne.

Główny problem to jak sprawdzić, że dana właściwość wywodzi się z innej. Nie ma czegoś takiego jak PropertyInfo.BaseProperty. Sprawdzenia dokonujemy dla metody ustawiającej i pobierającej właściwości. Przy czym nie wystarczy zwykłe porównanie z metodą bazową. Może się zdarzyć jeśli użyjemy właściwości automatycznych, że w hierarchii klas A->B->C, gdzie A ma właściwość abstrakcyjną, w B i C ją nadpisujemy to metody B i C jako bazową wskazują na A.

Sprawdzenie czy dwie właściwości są takie same też nie jest proste. Ze względu na możliwe różne ReflectedType porównanie wprost nie działa tak jak powinno. Porównujemy więc nazwę, typ właściwości i jej parametry indeksowe. Porównania parametrów indeksowych znowu nie możemy wprost, gdyż mogą się one różnić wartością MetadataToken. Porównujemy więc ich typy. Ponieważ do porównania nie bierzemy wszystkiego tylko kilka elementów nie jestem pewny na 100%, czy są jakieś przypadki które będą zakwalifikowane jako false positive albo false nagative.

private static IEnumerable<methodinfo> GetBaseDefinitions(this MethodInfo a_mi,
    bool a_with_this = false)
{
    if (a_with_this)
        yield return a_mi;

    MethodInfo t = a_mi;

    while ((t.GetBaseDefinition() != null) && (t.GetBaseDefinition() != t))
    {
        t = t.GetBaseDefinition();
        yield return t;
    }
}

private static bool ContainsAny<t>(this IEnumerable<t> a_enumerable, 
    IEnumerable<t> a_values)
{
    return a_enumerable.Intersect(a_values).Any();
}

private static bool IsAbstract(this PropertyInfo a_pi)
{
    if (a_pi.GetAccessors(true).Length == 0)
        return false;
    return a_pi.GetAccessors(true)[0].IsAbstract;
}

private static bool IsDerivedFrom(this PropertyInfo a_pi, PropertyInfo a_base,
    bool a_with_this = false)
{
    if (a_pi.Name != a_base.Name)
        return false;
    if (a_pi.PropertyType != a_base.PropertyType)
        return false;
    if (!a_pi.GetIndexParameters().Select(p => p.ParameterType).SequenceEqual(
        a_base.GetIndexParameters().Select(p => p.ParameterType)))
    {
        return false;
    }
    if (a_pi.DeclaringType == a_base.DeclaringType)
        return a_with_this;

    MethodInfo m1 = a_pi.GetGetMethod(true);
    MethodInfo m3 = a_base.GetGetMethod(true);

    if ((m1 != null) && (m3 != null))
    {
        if (m1.GetBaseDefinitions().ContainsAny(m3.GetBaseDefinitions(true)))
            return true;
    }
    else if ((m1 != null) || (m3 != null))
        return false;

    MethodInfo m2 = a_pi.GetSetMethod(true);
    MethodInfo m4 = a_base.GetSetMethod(true);

    if ((m2 != null) && (m4 != null))
    {
        if (m2.GetBaseDefinitions().ContainsAny(m4.GetBaseDefinitions(true)))
            return true;
    }
    else if ((m2 != null) || (m4 != null))
        return false;

    return false;
}

private static IEnumerable<type> GetBaseTypes(this Type a_type,
    bool a_with_this = false)
{
    if (a_with_this)
        yield return a_type;

    Type t = a_type;

    while (t.BaseType != null)
    {
        t = t.BaseType;
        yield return t;
    }
}

public static IEnumerable<propertyinfo> GetAllProperties(this Type type)
{
    List<propertyinfo> result = new List<propertyinfo>();

    foreach (var t in type.GetBaseTypes(true))
    {
        if (t == typeof(Object))
            break;
        if (t == typeof(ValueType))
            break;

        PropertyInfo[] type_props = t.GetProperties(
            BindingFlags.Public | BindingFlags.NonPublic |
            BindingFlags.Instance);

        foreach (var poss_prop in type_props.Reverse())
        {
            if (poss_prop.IsAbstract())
                continue;

            if (result.All(prop => !prop.IsDerivedFrom(poss_prop, true)))
                result.Add(poss_prop);
        }
    }

    result.Reverse();
    return result;
}

public static IEnumerable<fieldinfo> GetAllFields(this Type type,
    bool a_filter_autogenerated = true)
{
    foreach (var t in type.GetBaseTypes(true).Reverse())
    {
        if (t == typeof(Object))
            continue;
        if (t == typeof(ValueType))
            continue;

        FieldInfo[] fields = t.GetFields(BindingFlags.Instance | BindingFlags.Public |
            BindingFlags.NonPublic);

        foreach (var field in fields)
        {
            if (field.IsDefined(typeof(CompilerGeneratedAttribute), false) &&
                a_filter_autogenerated)
            {
                continue;
            }

            if (field.DeclaringType == t)
                yield return field;
        }
    }
}

Sprawdzenie czy właściwość wywodzi się z innej właściwości.

Metoda ta była mi potrzebna do wyciągnięcia z klasy wszystkich właściwości. W tym celu trzeba przejść przez wszystkie podtypu typu dla którego wyciągamy informację, dla każdego z nich trzeba pobrać wszystkie właściwości i wykluczyć powtarzające.

Niestety PropertyInfo niosą bardzo mało informacji w sobie i trzeba odwoływać się do ich metod getter- i setter-a.

W pierwszym etapie sprawdzamy, czy te właściwości nie są przypadkiem takie same. Niestety z uwagi na ReflectedType nie możemy tak po prostu je porównać.

Dalej analizujemy zależności pomiędzy ich metodami ustawiającymi i pobierającymi. Do analizy zależności wykorzystujemy MetohodInfo.GetBaseDefinition (są dwie przeładowane wersje, ta druga z parametrem jest dla nas bardziej użyteczna, gdyż zwraca niepubliczne metody). I tutaj znowu niestety. Jeśli automatyczna właściwość jest zdefiniowana w klasie bazowej jako wirtualna i nadpisana w dwóch kolejnych klasach, to obie nadpisane właściwości wskazują poprzez MetohodInfo.GetBaseDefinitions swoich metod na metody właściwości w klasie wirtualnej. I to zmusiło mnie z prostego sprawdzania, czy MetohodInfo.GetBaseDefinitions jednej jest równe GetBaseDefinitions drugiej do sprawdzenia czy nie mają one jakiejś wspólnej metody bazowej.

Kod:
private static IEnumerable GetBaseDefinitions(this MethodInfo a_mi,
    bool a_with_this = false)
{
    if (a_with_this)
        yield return a_mi;
        
    MethodInfo t = a_mi;
        
    while ((t.GetBaseDefinition() != null) && (t.GetBaseDefinition() != t))
    {
        t = t.GetBaseDefinition();
        yield return t;
    }
}

private static bool ContainsAny(this IEnumerable a_enumerable, 
    IEnumerable a_values)
{
    return a_enumerable.Intersect(a_values).Any();
}

private static bool IsDerivedFrom(this PropertyInfo a_pi, PropertyInfo a_base,
    bool a_with_this = false)
{
    if (a_pi.Name != a_base.Name)
        return false;
    if (a_pi.PropertyType != a_base.PropertyType)
        return false;
    if (!a_pi.GetIndexParameters().SequenceEqual(a_base.GetIndexParameters()))
        return false;
    if (a_pi.DeclaringType == a_base.DeclaringType)
        return a_with_this;
        
    MethodInfo m1 = a_pi.GetGetMethod(true);
    MethodInfo m3 = a_base.GetGetMethod(true);
        
    if ((m1 != null) && (m3 != null))
    {
        if (m1.GetBaseDefinitions().ContainsAny(m3.GetBaseDefinitions(true)))
            return true;
    }
    else if ((m1 != null) || (m3 != null))
        return false;
        
    MethodInfo m2 = a_pi.GetSetMethod(true);
    MethodInfo m4 = a_base.GetSetMethod(true);
        
    if ((m2 != null) && (m4 != null))
    {
        if (m2.GetBaseDefinitions().ContainsAny(m4.GetBaseDefinitions(true)))
            return true;
    }
    else if ((m2 != null) || (m4 != null))
        return false;
        
    return false;
}

Ustalanie modyfikatora właściwości i metod

Metody i własciwości mogą być oznaczone słowami kluczowymi virtual, override, abstract. Jak wykorzystując refleksję dobrać się do tych danych ?

Nie jest to niestety takie proste. Dla metod mamy dostępne właściwości MethodInfo.IsVirtual, MethodInfo.IsAbstract.

Informacje o właściwościach są bardzo ubogie. Całe wnioskowanie musimy przeprowadzić w oparciu o metody PropertyInfo.GetGetMethod i PropertyInfo.GetSetMethod do których dostęp możemy też uzyskać poprzez PropertyInfo.GetAccessors. Wszystkie trzy metody mają przeładowaną wersję za pomocą których dobieramy się do metod niepublicznych. Tak więc analiza właściwości sprowadza się do analizy jednej z jej metod: getter-a albo setter-a. Nie trzeba analizować obu gdyż modyfikatory typu virtual, override, abstract są nakładane na właściwość.

Do ustalenia metody bazowej służy MethodInfo.GetBaseDefinition. W przypadku braku takowej funkcja nie zwraca null tylko samą siebie.

Do ustalenia typu w którym metoda została zdefiniowana służy MethidInfo.DeclaringType. Jest jeszcze MethidInfo.ReflectedType, który zwraca typ z którego dana metoda została pobrana. Np. jeśli jest to metoda publiczna zdefiniowana w klasie bazowej A, a pobrana z klasy podrzędnej B to te dwie wartości będą różne. Co więcej porównanie metod da inne wyniki.

Całość ogólnie tworzy mętny obraz. Sam nie wiem czy podane niżej rozwiązanie obejmuje wszystkie przypadki i czy nie ma jakiś wyjątków. Cała refleksja w C# jest potężna, ale jak w nią wnikamy pojawia się mnóstwo problemów. Brakuje podanych niżej metod. Metod na porównanie czy elementy są takie same, niezależnie od różnicy ReflectedType, brakuje rozbudowanej PropertyInfo, tak by nie trzeba było analizować jej metody getter-a i setter-a. Bardzo trudno ustalić, że metoda interfejsu została zadeklarowana bez użycia virtual. Bardzo trudno ustalić czy metoda została zadeklarowana przy pomocy new. Nie da się w żaden sposób powiązać właściwości automatycznej z polem prywatnym. I pewnie jeszcze wiele innych.

Moja wersja metod, która na obecny czas działają:
/// 
/// With virtual keyword. Also includes interface implementations even without 
/// virtual keyword.
/// 
/// /// 
public static bool IsVirtual(this MethodInfo a_mi)
{
    return a_mi.IsVirtual && !a_mi.IsAbstract && !a_mi.IsOverriden();
}

/// 
/// With abstract keyword.
/// 
/// /// 
public static bool IsAbstract(this MethodInfo a_mi)
{
    return a_mi.IsAbstract;
}

/// 
/// With override keyword.
/// 
/// /// 
public static bool IsOverriden(this MethodInfo a_mi)
{
    return a_mi.DeclaringType != a_mi.GetBaseDefinition().DeclaringType;

}

/// 
/// With virtual keyword. Also includes interface implementations 
/// even without virtual keyword.
/// 
/// /// 
public static bool IsVirtual(this PropertyInfo a_pi)
{
    if (a_pi.GetAccessors(true).Length == 0)
        return false;

    return a_pi.GetAccessors(true)[0].IsVirtual();
}

/// 
/// With abstract keyword.
/// 
/// /// 
public static bool IsAbstract(this PropertyInfo a_pi)
{
    if (a_pi.GetAccessors(true).Length == 0)
        return false;
    return a_pi.GetAccessors(true)[0].IsAbstract();
}

/// 
/// With override keyword.
/// 
/// /// 
public static bool IsOverriden(this PropertyInfo a_pi)
{
    if (a_pi.GetAccessors(true).Length == 0)
        return false;

    return a_pi.GetAccessors(true)[0].IsOverriden();
}

2011-08-02

Mapowanie tekstury na sferę

Rysunek:


  
    
      
    
    
      
    
    
      
    
    
      
    
  
  
  
    
      
        image/svg+xml
        
        
      
    
  
  
    
    
      
      
    
    
      
      
    
    
      
      
    
    
      
      
    
    
    
    
      
      
    
    
      
      
    
    
    
    P
    O
    U
    F
    R
    A
    B
    φ
    ϴ
    
    
  
Wektory $\mathbf{\hat U}$ (góra), $\mathbf{\hat F}$ (wprzód), $\mathbf{\hat R}$ (prawo) orientując sferę w przestrzeni. Zakładamy, że wektory te są znormalizowane. Mapowanie przeprowadzamy dla punktu P znajdującego się na sferze.

Do wyznaczenia jest para (u,v), gdzie u zależy od kąta $\mathbf{\theta}$, v zależy od kąta $\mathbf{\varphi}$, $\mathbf{u\in \left < 0,1 \right >}$, $\mathbf{v\in \left < 0,1 \right >}$, $\mathbf{0\leqslant u+v\leqslant 1}$.

Poniżej podane wzory znacząco się upraszczają jeśli wektory $\mathbf{\hat U}$, $\mathbf{\hat F}$, $\mathbf{\hat R}$ są wersorami układu współrzędnych, a środek sfery to środek układu współrzędnych. Wtedy korzystamy z klasycznych wzorów na konwersję z układu kartezjańskiego na sferyczny.

Normalna (wektor znormalizowany) w punkcie P ma wartość:

$\hat N = \widehat{OP}$
Wyznaczenie wartości $\mathbf{\varphi}$:

$\cos \varphi = \hat U \cdot \hat N$

$\varphi = arccos \left ( \hat U \cdot \hat N \right )$
Wyznaczenie wartości $\mathbf{\theta}$:

$\vec A = \hat U \times \hat N$

$\hat B = \widehat{\vec A \times \hat U}$

$\cos \theta = \hat N \cdot \hat B$

$\theta = arccos \left ( \hat N \cdot \hat B \right )$
Mając wyznaczone oba kąty trzeba je zmapować na (u,v).

Chcemy by v zależał kąta $\mathbf{\varphi}$. Zauważmy że wartości funkcji arccos zmieniają się w zakresie $\mathbf{\left < 0, \pi \right >}$. Dla ujemnych wartości $\mathbf{\varphi}$ wartość arccos nie zmienia się co do znaku. Tak więc:

$v = \displaystyle\frac{\varphi}{\pi}$
W przypadku mapowania v znak kąta ma znaczenie. Chcemy by v zmieniało się w granicach od zera do jeden dla $\mathbf{\theta}$ zmieniającej się od 0 do $2\pi$. Ponieważ arccos zmienia się w zakresie $\left < 0, \pi \right >$ musimy rozpoznać czy punkt P leży z tyłu, czy z przodu sfery. W tym celu badamy znak $\mathbf{\cos \alpha = \hat F \cdot \hat B}$. Jest on równy kosinusowi kąta między tymi dwoma wektorami. Jeśli $\cos \alpha > 0$ to $\displaystyle-\frac{\pi}{2}<\alpha<\displaystyle\frac{\pi}{2}$. Czyli punkt P znajduje się z tyłu sfery. Dla punktów z tyłu sfery współrzędną u liczymy tak:

$v = \displaystyle\frac{\theta}{2\pi}$
Z przodu:

$v = 1 - \displaystyle\frac{\theta}{2\pi}$
Przykład implementacji w C#:
public override Vector2 GetUV(Intersection a_intersection)
{
    Vector3 normal = (a_intersection.Pos - Pos).Normalized;

    Vector3 B = Vector3.CrossProduct(
        Vector3.CrossProduct(Up, normal), Up).Normalized;

    double a = Forward * B;

    if (a < 0)
    {
        return new Vector2(
            1 - (Math.Acos(Right * B) / (2 * Constants.PI)),
            Math.Acos(Up * normal) / Constants.PI);
    }
    else
    {
        return new Vector2(
            Math.Acos(Right * B) / (2 * Constants.PI),
            Math.Acos(Up * normal) / Constants.PI);
    }
}
W przypadku gdy punkt na sferze do lokalnego układu współrzędnych obliczenia znacząco się upraszczają:
public virtual Vector3 GetUVW(Intersection a_intersection)
{
    return WorldToLocal * a_intersection.Pos;
}

public override Vector2 GetUV(Intersection a_intersection)
{
    Vector3 uvw = a_intersection.UVW;

    double v = Math.Acos(uvw.Y / uvw.Length) / Constants.PI;
    double u1 = Math.Atan2(uvw.Z, uvw.X) / (2 * Constants.PI);

    if (uvw.Z > 0)
        return new Vector2(u1, v);
    else
        return new Vector2(1 + u1, v);
}

Dodanie skalowanie obiektów

Dodanie skalowanie umożliwi nam skalowanie obiektów, dzięki czemu będziemy w stanie przeskalować element lub grupę elementów sceny bez modyfikacji ich współrzędnych i parametrów. Dzięki skalowaniu będziemy też mogli zrobić z sfery elipsoidę.

Dodanie skalowania wymusi zmiany na większości już napisanego kodu. Większość parametrów intersekcji i samego promienia przed wykorzystaniem w kodzie obiektów sceny będzie musiała zostać sprowadzona do lokalnego układu współrzędnych.

Dzięki temu kod w wielu miejscach będzie można uprościć gdyż obiekt w lokalnym układzie współrzędnych będzie zorientowany wzdłuż wersorów układu współrzędnych, a jego środek będzie w środku układu współrzędnych.

Z dużym prawdopodobieństwem zmiany wprowadzone przy okazji dodawania skalowania powinny także pociągnąć pochylanie.

Współczynniki skalowania mogą być także ujemne.

Orientacja obiektu w przestrzeni - skalowanie

Po dodaniu skalowania do informacji orientującej obiekt w przestrzeni kod tworzący macierze transformacji świat-obiekt zmienia się z :
private void UpdateTransformationMatrices()
{
    m_local_to_world = Matrix4.CreateTranslation(Pos) *
        new Matrix4(Right, Up, Forward);
    m_world_to_local = m_local_to_world.Inverted;
}
na:
private void UpdateTransformationMatrices()
{
    m_local_to_world = Matrix4.CreateTranslation(Pos) * 
        new Matrix4(Right, Up, Forward) * Matrix4.CreateScale(Scale);
    m_world_to_local = m_local_to_world.Inverted;
}

2011-07-27

Intersekcja promienia z trójkątem

Równanie promienia w postaci parametrycznej:

$P = S + \vec{d}t$
Współrzędne punkt trójkąta 3D, mając jego współrzędne barycentryczne liczymy:

$P=P_0+u(P_1-P_0)+v(P_2-P_0)$
Z porównania:

$S + \vec{d}t = P_0+u(P_1-P_0)+v(P_2-P_0)$
Porządkując względem niewiadomych t, u, v:

$S + \vec{d}t = P_0 + uP_1 - uP_0 + vP_2 - vP_0$

$t\vec{d} - u(P_1-P_0) - v(P_2-P_0)=P_0-S$

$-t\vec{d} + u(P_1-P_0) + v(P_2-P_0)=S-P_0$

$\left [ -\vec{d}, P_1-P_0, P_2-P_0 \right ] \begin{bmatrix}t \\ u \\ v \end{bmatrix}=S-P_0$
Oznaczamy $E_1=P_1-P_0$, $E_2=P_2-P_0$:

$\left [ -\vec{d}, E_1, E_2 \right ] \begin{bmatrix}t \\ u \\ v \end{bmatrix}=S-P_0$
Korzystając ze wzorów Cramera rozwiązanie ma postać:

$\begin{bmatrix}t \\ u \\ v \end{bmatrix}=
\displaystyle\frac{1}{\left | -\vec{d}, E_1, E_2 \right |}
\begin{bmatrix}
\left | T, E_1, E_2 \right | \\
\left | -\vec{d}, T, E_2 \right | \\
\left | -\vec{d}, E_1, T \right |
\end{bmatrix}$
Gdzie wyrażenia typu $\left | T, E_1, E_2 \right |$ oznaczają wyznacznik macierzy w której trzy wiersze definiują wektory $T, E_1, E_2$.

Korzystając z właściwości iloczynu mieszanego i iloczynu wektorowego możemy zapisać:

$\begin{bmatrix}t \\ u \\ v \end{bmatrix}=\displaystyle\frac{1}{E_1 \cdot (\vec{d} \times E_2)}\begin{bmatrix}
E_2 \cdot (T \times E_1) \\
T \cdot (\vec{d} \times E_2) \\
\vec{d} \cdot (T \times E_1)
\end{bmatrix}$
Oznaczamy $P=\vec{d}\times E_2$ i $Q=T\times E_1$:

$\begin{bmatrix}t \\ u \\ v \end{bmatrix}=\displaystyle\frac{1}{E_1 \cdot P}\begin{bmatrix}
E_2 \cdot Q \\
T \cdot P \\
\vec{d} \cdot Q
\end{bmatrix}$
Sprawdzając znak t upewniamy się czy przecięcie nie nastąpiło przed punktem startu promienia S.

Za pomocą u i v sprawdzamy czy punkt przecięcia nastąpił we wnętrzu trójkąta. Otrzymane niejako przy okazji współrzędne barycentryczne (u,v) przydają się podczas np. teksturowania, warto je sobie gdzieś zapamiętać, by nie musieć ich liczyć ponownie.

Jeśli $E_1 \cdot P=0$ układ nie ma rozwiązania. Promień biegnie równolegle do płaszczyzny w której zawiera się trójkąt. Wartość tego wyrażenia jest długością rzutu wektora kierunkowego prostej na normalną trójkąta. Badając znak $E_1 \cdot P=E_1$ ustalamy z której strony trójkąta uderzył w niego promień.

Przykład implementacji w C#:

public override Intersection GetIntersection(Ray a_ray)
{
    if (a_ray.PrevHit == this)
        return Scene.NoIntersection;

    if (OneSide)
    {
        if ((a_ray.Dir * Normal) > 0)
            return Scene.NoIntersection;
    }

    var E1 = V2 - V1;
    var E2 = V3 - V1;
       
    Vector3 P = Vector3.CrossProduct(a_ray.Dir, E2);
    double a = E1 * P;
            
    if (a == 0)
        return Scene.NoIntersection;
            
    Vector3 T = a_ray.Start - V1;
    double u = (P * T) / a;

    if (u < 0)
    {
        if (u < -INTERSECTION_ERROR) 
            return Scene.NoIntersection;
        u = 0;
    }

    if (u > 1)
    {
        if (u > 1 + INTERSECTION_ERROR)
            return Scene.NoIntersection;
        u = 1;
    }
            
    Vector3 Q = Vector3.CrossProduct(T, E1);
    double v = (a_ray.Dir * Q) / a;
            
    if (v < 0)
    {
        if (v < -INTERSECTION_ERROR) 
            return Scene.NoIntersection;
        v = 0;
    }
    if (u + v > 1.0)
        return Scene.NoIntersection;
            
    double dist = (Q * E2) / a;
            
    if (dist > 0)
    {
        bool backHit = a < 0;
            
        if (backHit && OneSide)
            return Scene.NoIntersection;
        else
        {
            return new Intersection()
            {
                SceneObject = this,
                SourceRay = a_ray,
                Dist = dist,
                Scene = Scene,
                BackHit = backHit,
                Pos = a_ray.HitPoint(dist),
                BaricentricUV = new Vector2(u, v),
                UV = BaricentricToUV(u, v)
            };
        }
    }
    else
        return Scene.NoIntersection;
}
Kontrola wartości (u,v) została tutaj zmodyfikowana o uwzględnienie błędów obliczeniowych. Dzięki czemu dwóm przystającym trójkątom nie zdarza się mieć prześwitów na przystającej krawędzi. Zastosowana w kodzie metoda na pewno zwiększa pole trójkąta, ale na tyle nieznacznie, że nie powinno mieć to znaczenia. Poza tym wymiana kolejności testowania ich intersekcji z promieniem spowoduje wygenerowanie (nieznacznie) różnej sceny. Stała INTERSECTION_ERROR została ustalona przeze mnie w kodzie testowym. Na chwile obecną jej wartość to 1e-13. I jest sto razy większa niż wyszło to w kodzie testowym. Na kod testowy składa się scena złożona z niebieskiego kwadratu, zaraz za nimi jest komplet żółtych kwadratów dokładnie pokrywający przekątną. Całą sceną losowo kręcimy lub przemieszczamy. Za każdym razem bierzemy nową scenę by uniknąć błędów obliczeniowych. Położenie kwadratu względem kamery nie zmienia. Próbkujemy więc na sztywno piksele wokół przekątnej i sprawdzamy czy są niebieskie. I tak z kilka tysięcy razy dla pewności.

2011-07-22

Mapowanie tesktury na płaszczyznę

Musimy ustalić współrzędne tekstury (u,v), gdzie $\mathbf{0<=u<=1}$, $\mathbf{0<=v<=1}$. Taką współrzędną łatwo później przeskalować na konkretną bitmapę lub wykorzystać w wyliczaniu tekstury proceduralnej.

Aby zorientować teksturę na płaszczyźnie potrzebujemy punktu na płaszczyźnie, który będzie początkiem układu współrzędnych i dwóch wektorów orientujących wokół tego punktu dwuwymiarowy układ współrzędnych leżący na płaszczyźnie.

Do zorientowania obiektu w 3D potrzebujemy jego pozycji i dwa wektory orientujące jego układ współrzędnych. Trzeci wektor możemy sobie wyliczyć z iloczynu wektorowego. Z punktu (przesunięcia) i 3 wektorów budujemy macierze transformacji świat-obiekt i obiekt-świat.

Mnożąc punkt należący do płaszczyzny przez macierz świat-obiekt sprowadzamy do układu współrzędnych obiektu. Zakładamy, że taki punkt ciągle leży na płaszczyźnie.

Jego współrzędne X i Z wyznaczą nam współrzędne tekstury (u,v).

Ostatni problem do rozwiązania to nieskończoność płaszczyzny, który wymusza na nas przeskalowanie współrzędnych tekstury do pożądanego przez nas wizualnego efektu.

Orientacja obiektu w przestrzeni

Każdy silnik czy to raytracingu czy renderingu umieszczając obiekty w układzie współrzędnych świata wykorzystuje informację o położeniu obiektu i jego orientacji.

Orientację wyznaczają trzy wektory, które wskazują w górę obiektu, w jego bok i naprzód (z reguł jednostkowe i wzajemnie prostopadłe). W tych wektorów i pozycji budujemy dwie macierze świat-obiekt i obiekt-świat, jedna jest odwrotnością drugiej.

Mnożąc punkt świata przez macierz świat-obiekt sprowadzamy go układu współrzędnych obiektu. Z reguły w takim układzie środek układu współrzędnych to środek obiektu.

Tworzenie macierzy transformacji:
private void UpdateTransformationMatrices()
{
    m_local_to_world = Matrix4.CreateTranslation(Pos) *
        new Matrix4(Right, Up, Forward);
    m_world_to_local = m_local_to_world.Inverted;
}
Zwróćmy uwagę, że przy transformacji świat-obiekt najpierw przesuwamy punkt, a następnie go obracamy. Przy transformacji odwrotnej obiekt-świat jest odwrotnie.

Często jest tak, że nasz obiekt zawiera w sobie inne, a one znowu inne, a każdy ma swoje macierze transformacji. Ostateczna macierz transformacji powstaje przez wymnożenie wszystkich składowych macierzy.

Mapowanie tekstury - informacje ogólne

Mając punkt uderzenia promienia w obiekt musimy ustalić współrzędne tekstury (u,v), gdzie $\mathbf{0<=u<=1}$, $\mathbf{0<=v<=1}$. Taką współrzędną łatwo później przeskalować na konkretną bitmapę lub wykorzystać w wyliczaniu tekstury proceduralnej.

2011-07-20

Intersekcja promienia z płaszczyzną

Równanie promienia w postaci parametrycznej wychodzącego z punktu P w kierunku $\mathbf{\vec d}$

$x = P + \vec{d}t$
Parametr t możemy potraktować jako czas. Jest ono identyczne z równaniem parametrycznym linii, nadajemy tylko innego znaczenia parametrom.

Równanie parametryczne płaszczyzny w postaci wektorowej:

$\vec n \cdot x + D = 0$
Z porównania:

$\vec n(P + \vec{d}t) + D = 0$
Stąd:

$t=\displaystyle-\frac{D + \vec{n}P}{\vec{n} \cdot \vec d}$
Jeśli $\mathbf{t<0}$ to przecięcie nastąpiło po przed miejscem wyjścia promienia. Traktując promień jako prostą takie t ma sens. Jeśli $\mathbf{t>=0}$ punkt przecięcia promienia z prostą istnieje.

Warunkiem istnienia punktu przecięcia jest to by płaszczyzn była nierównoległa do promienia, czyli $\mathbf{\vec n \cdot \vec d} \neq 0$.

Punkt przecięcia wyznaczamy podstawiając t do równania prostej.

Aby ustalić z której strony płaszczyzny nastąpiło uderzenie badamy znak wyrażenia $\mathbf{\vec{n} \cdot \vec d}$. Jeśli $\mathbf{\vec{n} \cdot \vec d}<0$ uderzenie nastąpiło z przodu płaszczyzny (z kierunku który wskazuje normalna płaszczyzny). Jeśli $\mathbf{\vec{n} \cdot \vec d}>0$ uderzenie nastąpiło z tyłu. Jeśli $\mathbf{\vec{n} \cdot \vec d}=0$ promień jest równoległy do płaszczyzny.

Przykład implementacji w C#:
public override Intersection GetIntersection(Ray a_ray)
{
    double denom = Normal * a_ray.Dir;

    if (OneSide)
    {
        if (denom > 0)
            return Scene.NoIntersection;
    }

    if (a_ray.PrevHit == this)
        return Scene.NoIntersection;

    if (!denom.IsZero())
    {
        double dist = ((Normal * a_ray.Start) - Distance) / (-denom);

        if (dist <= 0)
            return Scene.NoIntersection;
        else
        {
            bool backHit = (denom > 0);

            if (backHit && OneSide)
                return Scene.NoIntersection;
            else
            {
                return new Intersection()
                {
                    SceneObject = this,
                    SourceRay = a_ray,
                    Dist = dist,
                    Scene = Scene,
                    BackHit = backHit,
                    Pos = a_ray.HitPoint(dist)
                };
            }
        }
    }
    else
        return Scene.NoIntersection;
}

2011-07-19

Odległość punktu od płaszczyzny

Mamy dowolny punkt płaszczyzny $\mathbf{x_0}$, jej normalną $\mathbf{\vec n}$. Należy znaleźć odległość punktu $\mathbf{x}$ od płaszczyzny. Jest to odległość minimalna, czyli odległość pomiędzy punktem $\mathbf{x}$, a takim punktem płaszczyzny, że wektor wyznaczony przez te dwa punkty jest równoległy do normalnej płaszczyzny.

Odległość ta jest niczym innym jak długością rzutu wektora $\mathbf{\vec v = x_0 - x}$ na wektor normalny. Czyli:

$d=|(\vec v \cdot \vec n) \hat n|$
$d=|((x-x_0) \cdot \vec n) \hat n|$
Jeśli wektor $\mathbf{\vec n}$ jest znormalizowany to:

$d=\vec v \cdot \vec n$
$d=(x-x_0) \cdot \vec n$
W szczególności jeśli punkt $\mathbf{x=(0,0,0)}$ to otrzymujemy odległość płaszczyzny od początku układu współrzędnych ze znakiem minus.

Znak $\mathbf{d}$ wskazuje nam po której stronie płaszczyzny leży punkt $\mathbf{x_0}$.

Równanie płaszczyzny

Równanie płaszczyzny ma postać:

$\vec N \cdot (X - X_0)$
, gdzie $\mathbf{\vec N=[A,B,C]}$ to normalna płaszczyzny, $\mathbf{X=(x,y,z)}$ to dowolny punkt płaszczyzny, $\mathbf{X_0=[x_0,y_0,z_0]}$ to znany punkt płaszczyzny. Iloczyn skalarny dwóch wektorów wzajemnie do siebie prostopadłych równy jest zeru. Równianie to więc spełniają takie $\mathbf{x}$, które leżą na płaszczyźnie prostopadłej do $\mathbf{\vec N}$ i przechodzącej przez punkt $\mathbf{X_0}$.

Podstawiając współrzędne:

$\vec N \cdot (X - X_0)$
$[A,B,C] \cdot [x-x_0,y-y_0,z-z_0]=0$
$Ax+By+Cz+(-Ax_0-By_0-Cz_0)=0$
$Ax+By+Cz+D=0$
Jest to algebraiczna postać płaszczyzny. D to dowolna stała. Możmy powiedzieć, że -D to rzut punktu $\mathbf{X_0}$ na wektor $\mathbf{\vec N}$. Jest to minimalna odległość płaszczyzny od początku układu współrzędnych. $A,B,C,D mogą przyjmować dowolne wartości, tak by $\mathbf{|\vec N|} \neq 0$.

Równanie płaszczyzny możemy też zapisać jako:

$\vec N \cdot (X - X_0)$
$\vec N \cdot X - \vec N \cdot X_0)$
$\vec N \cdot X + D$
Bardzo często musimy określić płaszczyznę mając jej trzy punkty. Punkty powinny być niewspółliniowe. Budujemy z nich dwa wektory równoległe do płaszczyzny, ich iloczyn wektorowy to normalna.

Mając dowolny punkt płaszczyzny $\mathbf{p}$ minimalna odległość do początku układu współrzędnych to rzut wektora $\mathbf{[p_x-0,p_y-0,p_z-0]}$ na wektor $\mathbf{\vec N}$. Pod warunkiem, że wektor normalny jest znormalizowany jest to ich iloczyn skalarny.

Prawo załamania

Mając dany wektor promienia padający na granice dwóch ośrodków $\mathbf{\vec i}$, współczynniki załamania obu ośrodków $\mathbf{n_1}$ i $\mathbf{n_2}$, normalną w punkcie padania $\mathbf{\vec N}$ wyznaczymy wektor załamany $\mathbf{\vec t}$.

Zakładmy, że wektor normalny $\mathbf{\vec N}$, padający $\mathbf{\vec i}$ jak i załamany $\mathbf{\vec t}$ są znormalizowane.


  
    
      
    
    
      
    
    
      
    
    
      
    
    
      
    
  
  
  
    
      
        image/svg+xml
        
        
      
    
  
  
    
    
    
    
    
    
      
      
    
    
      
      
    
    
      
      
    
    
      
      
    
    
      
      
    
    
      
      
    
    
    α
    β
    
    
    
      
        t
        1
      
      
    
    
      t
      
    
    
      i
      1
      
    
    
      
        N
      
      
    
    
      
        t
        2
      
      
    
    
      2
      
        i
        
      
    
    
      i
      
    
    
      n
      1
    
    
      n
      2
    
  


Ponieważ wektory $\mathbf{\vec i}$ i $\mathbf{\vec t}$ są znormalizowane to:

$\left | \vec{i_2} \right | = \sin \alpha \qquad\qquad
\left | \vec{i_1} \right | = \cos \alpha \qquad\qquad
\left | \vec{t_2} \right | = \sin \beta \qquad\qquad
\left | \vec{t_1} \right | = \cos \beta$
Z rzutów wektorów $\mathbf{\vec i}$ i $\mathbf{\vec t}$ na $\mathbf{\vec N}$:

$\vec i \vec N = -\cos \alpha \qquad\qquad
\vec t \vec N = -\cos \beta$
Znając długości wektorów $\mathbf{\left | \vec{i_1} \right |}$ i $\mathbf{\left | \vec{t_1} \right |}$ i ich zwrot możemy zapisać:

$\vec{i_1} = -\vec N \cos \alpha \qquad\qquad
\vec{t_1} = -\vec N \cos \beta$
Z prawa Snelliusa:

$n_1 \sin \alpha = n_2 \sin \beta$

$n_1 \left | \vec{i_2} \right |= n_2 \left | \vec{t_2} \right |$
Ponieważ wektory $\mathbf{\vec{i_2}}$ i $\mathbf{\vec{t_2}}$ mają zgodny zwrot i kierunek:

$n_1 \vec{i_2} = n_2 \vec{t_2}$

$n_1 ( \vec i - \vec i_1 ) = n_2 \vec{t_2}$

$n_1 ( \vec i + \vec N \cos \alpha ) = n_2 \vec{t_2}$

$\vec{t_2} = \displaystyle\frac{n_1}{n_2}(\vec i + \vec N \cos \alpha)$
Teraz wyznaczmy $\mathbf{\vec{t_1}}$:

$\begin{split}\vec{t_1} &= -\vec N \cos \beta \\
&= -\vec N \left | \vec t_1 \right | \\
&= -\vec N \sqrt{1-\left | \vec t_2 \right |^2} \\
&= -\vec N \sqrt{1-\sin^2 \beta} \\
&= -\vec N \sqrt{1- \displaystyle\frac{n_1}{n_2}^2 \sin^2 \alpha} \\
&= -\vec N \sqrt{1- \displaystyle\frac{n_1}{n_2}^2 ( 1 - \cos^2 \alpha )}
\end{split}$
Mając $\mathbf{\vec{t_1}}$ i $\mathbf{\vec{t_2}}$:

$\begin{split}\vec t &= \vec{t_1} + \vec{t_2} \\
&= \displaystyle\frac{n_1}{n_2}(\vec i + \vec N \cos \alpha) - \vec N \sqrt{1- \displaystyle\frac{n_1}{n_2}^2 ( 1 - \cos^2 \alpha )} \\
&= \displaystyle\frac{n_1}{n_2}\vec i + \vec N ( \displaystyle\frac{n_1}{n_2} \cos \alpha - \sqrt{1- \displaystyle\frac{n_1}{n_2}^2 ( 1 - \cos^2 \alpha )})
\end{split}$
Otrzymany wektor $\mathbf{\vec{t}}$ jest znormalizowany.

Promień padający po w kierunku normalnej nie zostanie załamany.

Jeśli $\mathbf{n_1 < n_2}$ wraz ze zwiększaniem kąta $\mathbf{\alpha}$ kąt $\mathbf{\beta}$ będzie się zwiekszał, ale wolniej. Jak $\mathbf{\alpha}$ będzie dążył do $90^\circ$ to $\mathbf{\beta}$ będzie dążył do $arcsin\frac{n_1}{n_2}$, kąta mniejszego od $90^\circ$. Jeśli $\mathbf{n_2 < n_1}$ wraz ze zwiększaniem kąta $\mathbf{\alpha}$ kąt $\mathbf{\beta}$ będzie się zwiększał, ale szybciej. Kąt $\mathbf{\beta}$ zmienia się zgodnie z $arcsin(\frac{n_1}{n_2}\sin\alpha$). Dla pewnego granicznego kąta $\mathbf{\alpha}$ wartość z której liczymy arcsin będzie większa niż jeden. Nastąpi to dla kąta granicznego $arcsin\frac{n_2}{n_1}$. Po przekroczeniu tego kąta promień załamany zanika i następuje całkowite wewnętrzne odbicie.

W naszym równaniu na wektor promienia załamanego zjawisko całkowitego wewnętrznego odbicia następuje jeśli:

$1- \displaystyle\frac{n_1}{n_2}^2 ( 1 - \cos^2 \alpha ) < 0$



W odniesieniu do raytracing'u należy jeszcze uwzględnić jak zmienia się energia promienia odbitego i załamanego. Po pierwsze zakładamy, że istnieje współczynnik pochłaniania energii promienia odbitego i załamanego. Samo podział energii (a także polaryzacji) dla dielektryków definiuje równanie Fresnel'a. Sam promień przechodząc przez medium powinien tracić energię (Prawo Beer'a albo współczynniki pochłaniania liniowego i wykładniczego).

W niektórych mediach współczynnik załamania może się zmieniać w sposób nieliniowy w danym materiale i tak też będzie poruszać się promień światła w materiale.

Przykład implementacji w C#:
public Ray RefractedRay
{
    get
    {
        Debug.Assert(Exists);

        if (m_refracted_ray == null)
        {
            double n = RefractionIndiceN1 / RefractionIndiceN2;
            double cos_alpha = -SourceRay.Dir * ShadeNormal;
            double cos_beta_2 = 1 - n * n * (1 - cos_alpha * cos_alpha);

            if (cos_beta_2 > 0)
            {
                Scene.Statistics.RaysRefracted++;
                Scene.Statistics.RaysCreated++;

                m_refracted_ray = new Ray()
                {
                    Dir = n * SourceRay.Dir +
                        (n * cos_alpha - Math.Sqrt(cos_beta_2)) * ShadeNormal,
                    SourceIntersection = this,
                    Start = Pos,
                    Depth = SourceRay.Depth + 1,
                };
            }
            else
            {
                Scene.Statistics.TotalInternalReflections++;
            }
        }

        return m_refracted_ray;
    }
}

2011-07-18

Współrzędne barycentryczne

Ograniczymy się tylko do wersji 2D i skupimy się tylko na zastosowaniu w trójkącie.

  
    
      
    
  
  
  
    
      
        image/svg+xml
        
        
      
    
  
  
    
    P1
    P0
    P2
    v
    u
    
    
    
    P
    
  

Współrzędne barycentryczne to para (u,v). Można powiedzieć, że jest to punkt w układzie współrzędnych którego osie to boki trójkąta. Konwersja pomiędzy układem barycentrycznym, a ortogonalnym:

$P=P_0+\vec u + \vec v$
$P=P_0+u(P_1-P_0)+v(P_2-P_0)$
$P=P_0+uP_1-uP_0+vP_2-vP_0$
$P=P_0(1-v-u)+P_{1}u+P_{2}v$
Widzimy, że punkt $\mathbf{P}$ jest liniową kombinacją $\mathbf{P_0}$, $\mathbf{P_1}$ i $\mathbf{P_2}$. Tak naprawdę współrzędne barcentryczne to trójka liczb, przy czym ta trzecia to $\mathbf{1-u-v}$.

Współrzędne barycentryczne wyznaczają punkt wewnątrz trójkąta jeśli $\mathbf{0<=u<=1}$ i $\mathbf{0<=v<=1}$ i $\mathbf{0<=v+u<=1}$. Konwersja z układu ortgonalnego na układ barycentryczny:

  
    
      
    
  
  
  
    
      
        image/svg+xml
        
        
      
    
  
  
    
    P1
    P0
    P2
    P
    
    
  
Współrzędne barycentryczne wyznaczamy z proporcji pól trójkątów:

$u=\displaystyle\frac{\triangle(P, P0, P2)}{\triangle(P0, P1, P2)}$


$v=\displaystyle\frac{\triangle(P, P0, P1)}{\triangle(P0, P1, P2)}$

Pola trójkątów najlepiej policzyć z iloczynu wektorowego. Wzory te działają tylko dla punktu wewnątrz trójkąta.

Przeanalizujmy współrzędna u. Jeśli znajdzie się ona na boku trójkąta $\mathbf{\overline{P0P1}}$ pole trójkąta $\mathbf{\triangle(P, P0, P2)}$ jest równe zeru i tym samym współrzędna u. Wraz z przemieszaniem się punktu P w stronę punktu P1 pole trójkąta $\mathbf{\triangle(P, P0, P2)}$ dąży do pola trójkąta $\mathbf{\triangle(P0, P1, P2)}$, czyli współrzędna u dąży do jeden. Jeśli współrzędna wyjdzie poza trójkąt pole (powinno być, w naszym przypadku nie jest ujemne) lub większe od 1. To samo tyczy się pozostałych współrzędnych barycentrycznych.

Jeśli trójkąty $\mathbf{\triangle(P, P0, P2)}$ lub $\mathbf{\triangle(P, P0, P1)}$ będą mieć w sumie pole większe od pola trójkąta $\mathbf{\triangle(P0, P1, P2)}$ to niewątpliwie punkt P musi być poza trójkątem. Czyli punkt jest w trójkącie jeśli $\mathbf{0<=v+u<=1}$. Ten warunek musi być spełniony dla dowolnej pary współrzędnych barycentrycznych.

Współrzędne barycentryczne mają zastosowanie podczas teksturowania. Dla punktu P wewnątrz trójkąta liczymy jego współrzędne barycentryczne. Dalej biorąc nową trójkę współrzędnych tekstury liczymy współrzędną tekstury z współrzędnych barycentrycznych.

Reguła prawej dłoni

Reguła ta wyznacza nam zwrot iloczynu wektorowego. Jeśli kciuk wskaże pierwszą pierwszy wektor iloczynu. Pozostałe palce drugi wektor. To z wnętrza dłoni będzie wychodził wektor wynikowy. Dłoń jest rozwarta. Wszystkie palce są w jednej płaszczyźnie. Kciuk jest ustawiony prostopadle do pozostałych palców.

Dla układów prawoskrętnych reguła także wyznacza nam wizualny zwrot wektora. Dla układów lewoskrętnych musimy do tego stosować lewą rękę. Niezależnie od układu współrzędnych iloczyn wektorowy liczymy tak samo (jego wartość nie zmienia się).

Dodatkowo kierunek od kciuka do pozostałych palców wyznacza nam obrót wokół wektora wychodzącego z dłoni przy dodatnich kątach obrotu.

Tutaj został podany jeden ze sposobów jej definiowania. Można wykorzystać zaciśniętą pięść. Dłoń z palcami ustawionymi jak osie układu współrzędnych 3D. Można tą regułę tak zmodyfikować, że to kciuk wyznacza nam wektor wynikowy.

Skrętność układu współrzędnych

Układ współrzędnych może być lewo albo prawoskrętny. Własność ta definiuje jak wizualizujemy nasz układ współrzędnych. W którą stronę pokażemy zwrot iloczynu wektorowego. Z tym zaś wiąże się kierunek obrotu wektora wokół innego wektora i jak mierzymy kąt między wektorami ze znakiem.

Właściwość $\mathbf{X \times Y = Z}$ nie zależy od rodzaju układu współrzędnych.

W układzie lewoskrętnym jeśli X wskazuje w prawo, Y do góry, to Z wskazuje wgłąb sceny. W układzie prawoskrętnym Z wskazuje w stronę obserwatora sceny. Układ prawoskrętny to np. OpenGL, układ lewoskrętny to np. DirectX.

Jeśli chodzi o kąty. Dla układu prawoskrętnego bierzemy prawą rękę. Dla lewoskrętnego lewą. Kciuk wyznacza nam wektor wokół którego chcemy obrócić inny wektor. Palce złożone w pięść wyznaczają kierunek obrotu dla dodatnich wartości kąta.

We wszystkich przykładach posługuję się układem lewoskrętnym.

Iloczyn mieszany wektorów

Jest to liczba definiowana jako:

$\vec A \cdot (\vec B \times \vec C)$
Jest wartość bezwzględna stanowi objętość równoległościanu definiowanego przez te trzy wektory. Jeśli wszystkie trzy wektory leża na jednej płaszczyźnie wtedy iloczyn mieszany ma wartość zero.

Właściwości:

$\vec A \cdot (\vec B \times \vec C) = \vec B \cdot (\vec C \times \vec A) = \vec C \cdot (\vec A \times \vec B)$

$\vec A \cdot (\vec B \times \vec C) = \begin{vmatrix}a_1 & a_2 & a_3 \\ b_1 & b_2 & b_3 \\ c_1 & c_2 & c_3 \end{vmatrix}$