Zufällige 3D Vektoren innerhalb einer Kugel

hustbaer

@pumuckl:
Ich glaube das nennt man Skybox

SeppJ

314159265358979 schrieb:

Ich habe dein Verfahren, SeppJ, leider nicht verstanden, denn wir haben Integralrechnung in der Schule auch noch nicht gemacht. Es wäre schön, wenn mir das jemand programmieren könnte.

Ich habe im Moment gerade nicht die Mittel grossartig zu coden und an meinen schon fertigen Code fuer das Kugelproblem komme ich von hier gerade auch nicht dran. Wenn es dir morgen reicht, dann kann ich es dann machen, falls ich dran denke.

314159265358979

Alles klar. Alleine, dass du mir deinen 10000ten Post widmest, ist mir eine Ehre

krümelkacker

Das einfachste ist natürlich X/Y/Z aus einem die Kugel umschließenden Würfel gleichverteilt zu ziehen und das solange zu wiederholen, bis der Radius passt.

Sonst könnte man das vielleicht noch so machen:
Richtung und Radius unabhängig ziehen.

Gleichverteilte Richtung:
X/Y/Z/ Gauß-verteilt ziehen und normalisieren

Radius ziehen:
Die Wahrscheinlichkeitsdichte für den Radius muss mit dem Radius quadratisch zunehmen, also:
P(r) \propto r^2
da die Kugeloberfläche ja mit dem Radius auch quadratisch zunimmt. Sonst bekommt man hinterher keine Gleichverteilung.

*Rewind*

Das grundsätzliche Problem bei dieser Aufgabe ist, dass der Gültigkeitsbereich (und damit auch das Maximum) jeder Dimension sich in Abhängigkeit von den Werten der anderen zwei Dimensionen ergibt.
Ich habe zwar nichts von dem verstanden, was SeppJ erklärt hat, aber aus dem genannten Grund sehe ich keine andere Möglichkeit als den Würfel.

krümelkacker schrieb:

Radius ziehen:
Die Wahrscheinlichkeitsdichte für den Radius muss mit dem Radius quadratisch zunehmen, also:
P(r) \propto r^2
da die Kugeloberfläche ja mit dem Radius auch quadratisch zunimmt. Sonst bekommt man hinterher keine Gleichverteilung.

Wie bastelt man das dann im Endeffekt? Wenn Radius zwischen r1 und r2 liegen soll: Kann man einfach gleichverteilt zwischen r1*r2 und r2*r2 ziehen und daraus die Wurzel als Radius nehmen? Müsste doch eine quadratische Verteilung ergeben.

hustbaer

[Rewind] schrieb:

Das grundsätzliche Problem bei dieser Aufgabe ist, dass der Gültigkeitsbereich (und damit auch das Maximum) jeder Dimension sich in Abhängigkeit von den Werten der anderen zwei Dimensionen ergibt.

Nur wenn du mit X/Y/Z Koordinaten rechnest.
Wenn du dagegen mit Kugelkoordinaten (quasi Längengrad/Breitengrad/Radius) rechnest, dann sind alle drei Werte unabhängig.

Einer der Winkel hat immer den Bereich [0°...180°), der andere immer [0...360°), und der Radius eben immer [0...Radius].

Wenn du dann noch wie krümelkacker & SeppJ beschreiben haben die Verteilung anpasst, dann bekommst du auch Punkte raus, die im X/Y/Z Raum gleichverteilt sind.

fdfdg schrieb:

krümelkacker schrieb:

Radius ziehen:
Die Wahrscheinlichkeitsdichte für den Radius muss mit dem Radius quadratisch zunehmen, also:
P(r) \propto r^2
da die Kugeloberfläche ja mit dem Radius auch quadratisch zunimmt. Sonst bekommt man hinterher keine Gleichverteilung.

Wie bastelt man das dann im Endeffekt? Wenn Radius zwischen r1 und r2 liegen soll: Kann man einfach gleichverteilt zwischen r1*r2 und r2*r2 ziehen und daraus die Wurzel als Radius nehmen? Müsste doch eine quadratische Verteilung ergeben.

Stell eine Wahrscheinlichkeitsdichte Proportional zu r^2 auf und bestimme den Faktor mit der Bedingung, dass das Integral der Dichte von r1 bis r2 gleich 1 wird. Dann wendet man die Inversionsmethode an.

http://de.wikipedia.org/wiki/Inversionsmethode

µ schrieb:

fdfdg schrieb:

krümelkacker schrieb:

Radius ziehen:
Die Wahrscheinlichkeitsdichte für den Radius muss mit dem Radius quadratisch zunehmen, also:
P(r) \propto r^2
da die Kugeloberfläche ja mit dem Radius auch quadratisch zunimmt. Sonst bekommt man hinterher keine Gleichverteilung.

Wie bastelt man das dann im Endeffekt? Wenn Radius zwischen r1 und r2 liegen soll: Kann man einfach gleichverteilt zwischen r1*r2 und r2*r2 ziehen und daraus die Wurzel als Radius nehmen? Müsste doch eine quadratische Verteilung ergeben.

Stell eine Wahrscheinlichkeitsdichte Proportional zu r^2 auf und bestimme den Faktor mit der Bedingung, dass das Integral der Dichte von r1 bis r2 gleich 1 wird. Dann wendet man die Inversionsmethode an.

http://de.wikipedia.org/wiki/Inversionsmethode

Oder kurz: Ja, das stimmt so, fdfdg. Bis auf den Typo und die exakten Exponenten halt

Wen es interessiert, hab hier kurz was gebastelt: http://ideone.com/IsAj2

Werner Salomon

otze schrieb:

Für manche Dinge braucht es nicht C++11 sondern einfache Mathematik.

genau .. und es kommt öfter vor, als man denkt

Folgender Code sollte in etwa das sein, was PI braucht. Er sollte Positionen in einer Kugelschale zwischen zwei gegebenen Radien R1 und R2 erzeugen, die im Volumen annähernd gleich verteilt sind.

Im Gegensatz zu pumuckl und fdfdg bin ich aber der Meinung, dass man beim Radius die 3.Potenz nehmen muss, da der Radius im Volumen mit der 3.Potenz eingeht. Bei einer Kugelschale mit den Innenradius 1 und Außenradius 2 befinden sich 50% des Volumens jenseits eines Radius von ca. 1,651 bzw. = 3.Wurzel( 1 + (2^3 - 1^3)/2 ).

#include <boost/random/uniform_real.hpp>
#include <boost/random/linear_congruential.hpp> // minstd_rand

#include <iostream>
#include <cmath>

class Vector // nur für Demozwecke
{
public:
    typedef double value_type;
    Vector( double x, double y, double z )
        : m_x(x), m_y(y), m_z(z)
    {}
    friend Vector operator*( double factor, const Vector& v )
    {
        return Vector( v.m_x * factor, v.m_y * factor, v.m_z * factor );
    }
    friend double abs( const Vector& v )
    {
        return std::sqrt( v.m_x * v.m_x + v.m_y * v.m_y + v.m_z * v.m_z );
    }
    friend std::ostream& operator<<( std::ostream& out, const Vector& v )
    {
        return out << v.m_x << " " << v.m_y << " " << v.m_z;
    }
private:
    double m_x, m_y, m_z;
};

namespace
{
    const double PI = std::acos(-1.0);
} 

// -- Verteilung zur Generierung eines 'gleichverteilten' Einheitsvektors
template< typename V, typename T = typename V::value_type >
class Sphere
{
public:
    typedef T input_type;
    typedef V result_type;

    Sphere() 
        : m_sin_phi( -1., +1. )
        , m_beta( -PI, PI )
    {}
    template< typename G >
    result_type operator()( G& gen )
    {
        using std::cos;
        using std::sin;
        using std::sqrt;
        // siehe auch http://de.wikipedia.org/wiki/Kugelkalotte
        const T sin_phi = m_sin_phi( gen ); // sin( Breitengrad ) ist gleichverteilt, da H=2*PI*sin(phi)
        const T cos_phi = sqrt( 1 - sin_phi*sin_phi );
        const T beta = m_beta( gen ); // Längenwinkel
        return result_type( cos( beta )*cos_phi, sin( beta )*cos_phi, sin_phi );
    }
private:
    boost::uniform_real< T > m_sin_phi;
    boost::uniform_real< T > m_beta;
};

// --   Verteilung zur Generierung eines Radius im Intervall [x0, x1)
//      mit dem Ziel im Volumen einer Kugelschale eine 'Gleichverteilung' zu erreichen
template< typename T >
class Cubic
{
public:
    typedef T input_type;
    typedef T result_type;
    Cubic( result_type x0, result_type x1  )
        : m_distribution3( x0 * x0 *x0, x1 * x1 * x1 )
    {}
    template< typename G >
    result_type operator()( G& gen )
    {
        using std::pow;
        return pow( m_distribution3( gen ), result_type( 1/3. ) );
    }
private:
    boost::uniform_real< T > m_distribution3;
};

int main()
{
    using namespace std;
    const double R1 = 1.5;
    const double R2 = 4;
    boost::minstd_rand gen;
    Sphere< Vector > sphere;
    Cubic< double > cubic( R1, R2 );
    cout << fixed;
    cout.precision( 2 ); // 2 NK reichen hier
    for( int i=0; i<20; ++i )
    {
        Vector pos = cubic( gen ) * sphere( gen );
        cout << pos << "  " << abs(pos) << endl;
    }
    return 0;
}

Gruß
Werner

@Edit: in den Zeilen 42 und 58 den Typ 'Vector' gegen das Template-Argument V bzw. result_type ersetzt

SeppJ

@Werner: Bist du da sicher, dass die Verteilung deiner Längen stimmt? Dein Code ist mir zu verschlungen, um ihn jetzt noch zu entziffern. Du ziehst nur die dritte Wurzel, oder? Das würde nämlich nicht stimmen, da du dann nämlich entweder die Verteilung bei R1 als 0 hast oder die Verteilung verschoben (und somit die falsche Steigung) hast, je nachdem, was du da genau machst (sorry, wenn ich gerade keine 100 Zeilen Templates dekodieren möchte - habe heute Abend schon genug Denksport betrieben).

Ich habe jetzt stundenlang an der exakten Verteilung der Längen rumgetüftelt, ich kann nicht glauben, dass du so einen einfachen Code heraus bekommst, ich hingegen dieses Monstrum:

#include <iostream>
#include <vector>
#include <functional>
#include <numeric>
#include <algorithm>
#include <cstdlib>
#include <cmath>
#include <complex>
using namespace std;

int main ()
{
  vector<double> bins(1000,0.);
  const double r1 = 2;
  const double r2 = 5;
  double r1_3 = r1*r1*r1;
  double len = r2-r1;
  double cut = r1*r1 / (pow(r2,3)/3. - r1_3/3.);  // Wert der Verteilung bei r1 -> Offset
   double fac1 = cut * len;   // Der Offset muss noch mit der Länge normiert werden
  for(unsigned i = 0; i < 1000000; ++i)
    {
      double a = rand() *1. / RAND_MAX;
      if (a < fac1)              // Gleichverteilungsanteil der Verteilung
        a = (a/fac1)*len + r1;   // Dies ist noch einfach...
      else                       // ...aber jetzt kommt die x^2 Funktion die da drüber gelegt ist: OMG!
        {
          a = (a-fac1) / (1-fac1); // a ist wieder Einheitsverteilung
          // Herleitung: Mathematica, Zettel und Stift und viel Geduld
          double fac2 = len*len*(2*r1+r2)*a;
          complex<double> fac3 = pow(-2*r1_3+fac2+sqrt(complex<double>(fac2*(-4*r1_3+fac2),0.)),1./3.);
          complex<double> result = 0.5 * (
                                          pow(2,2./3.)*fac3
                                          +2*r1*
                                          (-1. + pow(2,1./3.)*r1 / fac3 )
                                  );
          a = result.real();   // Zusatzwissen: Imaginärteile heben sich genau weg
          a += r1;
        }
      ++bins[a*100];
    }
  double factor = 1./(accumulate(bins.begin(), bins.end(), 0.)*0.01);
  transform (bins.begin(), bins.end(), bins.begin(), bind1st(multiplies<double>(),factor));
  for (unsigned i = 0; i < bins.size(); ++i)
    {
      cout << i / 100. << '\t' << bins[i] << '\n';
    }
}

Dieser Code ist erfolgreich geprüft (beziehungweise er ist mein Prüfcode, einfach mal die Werte plotten, wer es nachvollziehen möchte) gegen die theoretische Wahrscheinlichkeitsverteilung von

P(R) = R^2/(R2^3/3 - R1^3/3) * Theta(R-R1) * Theta(R2-R)

Irgendwie muss man auch die komplexwertige Zwischenrechnung vereinfachen können, aber das würde ich selber gerne sehen und Mathematica schafft es nicht.

Die Längen kann man dann mit Punkten auf einer Einheitssphäre kombinieren, wie hier schon gezeigt wurde. Herleitung findet man hier:
http://mathworld.wolfram.com/SpherePointPicking.html
Das ist zum Glück ein ganzes Ende einfacher.

Ob sich die Rechnung nun lohnt sei dahingestellt. Bei der Wegwerfmethode kann man höchstens (4/3 r^3)/(2r)3 = 1/6 der Zufallszahlen benutzen, wenn man eine Vollkugel hat. Bei einer dünnen Schale sinkt der Anteil sehr schnell nahe 0. Und man hat natürlich auch milde Rechenkosten für die Abstandsberechnung.

Entschuldigung auch an alle Leser, dass mein Code ebenso unleserlich ist wie Werners, wenn auch aus weitaus weniger künstlerisch anspruchsvollen Gründen :p

P.S.: Herleitung ohne funktionierende LaTeX-Tags möchte ich niemandem antun. Im Prinzip ist es nichts kompliziertes, bloß eklig lang ohne CAS (und mit CAS auch ). Die Grundzüge wurden hier im Thread schon erklärt oder Verlinkt, das einzige was glaube ich noch nicht genannt wurde ist, wie man das macht, wenn die Inverse der kumulierten Wahrscheinlichkeitsverteilung bei 0 nicht 0 ist. Das ist das was hier im ersten if-Zweig behandelt wird, indem eine Gleichverteilung passender Höhe "abgeschröpft" und einfach (nach ein bisschen Skalieren) durchgeleitet wird. Der zweite if-Zweig ist dann der Teil der übrig bleibt und mit herkömmlichen Mitteln behandelt werden kann (auch wenn die herkömmlichen Mittel hier ganz schön komplexe Rechnung erfordern). Die Tatsache, dass die Imaginärteile bei negativer Wurzel sich gegenseitig weg heben habe ich zugegebenermaßen nur durch Ausprobieren getestet, aber konstruktionsbedingt muss sich das auch so ergeben, weil die Methodik an sich bewiesen ist und in der Rechnung selbst bis dahin keine Fehler sind (weil der Computer gerechnet hat ).

krümelkacker

SeppJ: Probier's doch mal so zu trenne, wie ich es vorgeschlagen habe. Wenn Du X/Y/Z Gauß-verteilt ziehst und dann normalisierst, hast Du schon auf der Kugeloverfläche eine Gleichverteilung. Dann fehlt nur noch der Abstand zum Ursprung. Ich habe den letzten Schritt nicht zu Ende gerechnet, aber so kompliziert dürfte das nicht sein.

*Rewind*

@Werner: wie lange hast du für den Code gebraucht?

otze

Bashar schrieb:

otze schrieb:

Für manche Dinge braucht e snicht C++11 sondern einfache Mathematik.

"einfache"? Ich hab grad keine Ahnung, wie das mathematisch gehen könnte.

ich meinte in der Tat die ziehen und Wegwerfmethode. Bei den typischen schnellen Zufallszahlengeneratoren ist die analytische Lösung wahrscheinlich langsamer als wegwerfen und neu ziehen. In dem Setup hier ist das auch gut möglich. Ist auch nicht weniger mathematisch, läuft dann halt unter Akzeptanzkriterium.

SeppJ

krümelkacker schrieb:

Ich habe den letzten Schritt nicht zu Ende gerechnet, aber so kompliziert dürfte das nicht sein.

Ja, das hatte ich auch gedacht:

SeppJ schrieb:

hier ist das bei der parabelfoermigen Verteilung der r-Koordinate aber noch einfach

Dankste .

Die Gleichverteilung auf der Kugelfläche ist nämlich der einfache Teil, auch wenn es zuerst nicht so aus sieht.

otze

SeppJ schrieb:

Die Gleichverteilung auf der Kugelfläche ist nämlich der einfache Teil, auch wenn es zuerst nicht so aus sieht.

Wenn wir problemlos für einen gezogenen Radius normalverteilte Punkte auf der Oberflöche erzeugne können, dann müssen wir doch nur den radius ordentlich ziehen.

Mein Ansatz wäre: wir teilen das Volumen der Kugel in unendlich viele Kugeloberflächen auf. dann wird ein Wert des Radius mit dem Verhätlnis gezogen, den seine (infinitesimal dicke) Kugeloberfläche am Integral über alle Oberflächen hat. Das ist praktischerweise das Volumen der Schicht

angenommen der maximale Radius ist 1 und a ist der innere Radius der Kugel, 0 < a < 1

also p(r) = 4* pi * r^2 /(4/3*pi - 4/3*pi *a^3) = 3*r^2/(1-a3), a <= r <= 1

darüber dann hübsch die Inversionmethode anwenden.

Sei q Normalverteilt in U(a,1)

q = \int_a^r' 3r^2/(1-a3) d r
<=> q = (r'^3 - a^3)/(1-a3)
<=> r = (q(1-a^3)+a3)^(1/3)

vermutlich habe ich mich irgendwo verrechnet.

Und dann 2 Winkel gleichverteilt ziehen, zusammen mit dem Radius müsste das einen gleichverteilten Vektor auf der Kugelschicht mit Radius r erzeugen.

Werner Salomon

SeppJ schrieb:

@Werner: Bist du da sicher, dass die Verteilung deiner Längen stimmt? Dein Code ist mir zu verschlungen, um ihn jetzt noch zu entziffern. Du ziehst nur die dritte Wurzel, oder? Das würde nämlich nicht stimmen, da du dann nämlich entweder die Verteilung bei R1 als 0 hast oder die Verteilung verschoben (und somit die falsche Steigung) hast, je nachdem, was du da genau machst

Ja ich ziehe die 3.Wurzel aus einer Verteilung von x0^3 bis x1^3 (s. Cubic<>) - das ist genau das, was otze vorschlägt

otze schrieb:

Mein Ansatz wäre: wir teilen das Volumen der Kugel in unendlich viele Kugeloberflächen auf. dann wird ein Wert des Radius mit dem Verhältnis gezogen, den seine (infinitesimal dicke) Kugeloberfläche am Integral über alle Oberflächen hat. Das ist praktischerweise das Volumen der Schicht

was dann kommt habe ich zwar nicht verstanden - insbesondere die Ausgangsgleichung, aber das Ergebnis

otze schrieb:

<=> r = (q*(1-a^3)+a3)^(1/3)

ist identisch mit dem, was die Verteilung class Cubic<> macht - siehe Zeile 74 und 80. Man ersetze a in otzes Ergebnis durch das Verhältnis x0/x1 aus Cubic<> und rechne es selber nach.

otze schrieb:

vermutlich habe ich mich irgendwo verrechnet.

nö - bis hierher anscheinend nicht.

otze schrieb:

Und dann 2 Winkel gleichverteilt ziehen, zusammen mit dem Radius müsste das einen gleichverteilten Vektor auf der Kugelschicht mit Radius r erzeugen.

warum das falsch ist, ist am Anfang von http://mathworld.wolfram.com/SpherePointPicking.html erklärt und durch eine Graphik demonstriert.

SeppJ schrieb:

Die Längen kann man dann mit Punkten auf einer Einheitssphäre kombinieren, wie hier schon gezeigt wurde. Herleitung findet man hier:
http://mathworld.wolfram.com/SpherePointPicking.html
Das ist zum Glück ein ganzes Ende einfacher.

Danke für den Link SeppJ,
Dort wird genau das beschrieben, was ich in der Klasse Sphere<> codiert habe, nur das die Bezeichner andere sind und die Intervalle der Winkel phasenverschoben sind - z.B. [0,2PI) statt [-PI,+PI). Aber das sollte das Ergebnis nicht verändern.

[Rewind] schrieb:

@Werner: wie lange hast du für den Code gebraucht?

weiß ich nicht mehr genau. Weniger als eine Stunde auf jeden Fall - vielleicht 45min, die meiste Zeit habe ich dafür gebraucht, aus der Formel für die Oberfläche der Kugelkalotte die Verteilung für den Einheitsvektor zu zaubern. Das ist das, was in den Zeilen 55 bis 57 steht.

Gruß
Werner

Werner Salomon

SeppJ schrieb:

Ich habe jetzt stundenlang an der exakten Verteilung der Längen rumgetüftelt, ich kann nicht glauben, dass du so einen einfachen Code heraus bekommst, ich hingegen dieses Monstrum:

Ich habe noch mal ein kleines Testprogramm geschrieben.
Es ist wirklich so, dass Dein 'Monstrum' exakt die selbe Verteilung liefert wie meine Verteilung class Cubic<>. Bei der identischen Zufallszahl liefern unsere Algorithmen zwar unterschiedliche Radien, aber bei näherer Betrachtung muss das auch so sein, aber die Verteilung ist die selbe.
Also Du hast Dich nicht verrechnet.

SeppJ schrieb:

Die Grundzüge wurden hier im Thread schon erklärt oder Verlinkt, ...

Das kann ich in keiner Weise nachvollziehen. Kannst Du da noch ein wenig konkreter werden, wie die Ausgangsformel oder Idee ist, die dann zu Deiner Lösung führt. Das würde mich interessieren.

Gruß
Werner

SeppJ

Werner Salomon schrieb:

SeppJ schrieb:

Die Grundzüge wurden hier im Thread schon erklärt oder Verlinkt, ...

Das kann ich in keiner Weise nachvollziehen. Kannst Du da noch ein wenig konkreter werden, wie die Ausgangsformel oder Idee ist, die dann zu Deiner Lösung führt. Das würde mich interessieren.

Im Nachhinein gesehen war's eher umständlich, auf otzes Trick hätte ich kommen müssen. Mein's ist im Prinzip die Inversionsmethode direkt mit dem Kopf voran durch die Betonwand gerechnet, anstatt die Tür zu nehmen .

Die Verteilung die ich haben möchte ist

P(R) = R^2/(R2^3/3 - R1^3/3) * Theta(R-R1) * Theta(R2-R)

Ich habe eine GLeichverteilung 0-1 die ich dazu umformen möchte. Wenn ich nun für die Inversionsmethode die kumulierte Wahrscheinlichkeitsverteilung zu berechnen, bekomme ich eine Funktion die bis R1 0 ist und dann wie x^3 bis 1 bei R2 steigt. Leider nicht direkt umkehrbar, weil mehrdeutig bei 0.

Und an dieser Stelle hätte ich es wie otze machen sollen . Was ich stattdessen gemacht habe: Das Problem, wie ich es sah, war, dass die Verteilung bei R1 nicht 0 ist. Wenn's einfach nur quadratisch von R1 bis R2 ginge und P(R1)=0, wäre das Problem einfach. Aber es ist im Prinzip ein "Block" von Gleichverteilung unter der quadratischen Verteilung. Daher habe ich gesagt, ich nehme aus meinen gleichverteilten Zufallszahlen einen Block eben dieser Größe (nämlich P(R1)*(R2-R1)) heraus und strecke diesen, so dass er diesem Block unter meiner Parabel entspricht. Das ist der erste if-Zweig bei mir.

Den Rest der gleichverteilten Zufallszahlen habe ich dann wieder gestreckt, so dass ich wieder eine Gleichvertelung habe, so dass ich die Inversionsmethode nach Schulbuch anwenden kann (Zeile 27). Nun hatte ich nur noch den parabelförmigen Teil von P(R) um den ich mich kümmern musste. Auf den Nullpjunkt verschoben und die Inversionsmethode direkt angewendet. Doch dummerweise haben sich in der Rechnung bis dahin derart viele Vorfaktoren angesammelt, dass das dann doch nicht mehr so einfach war, wie der Plan vorsah. Als kumulierte Wahrscheinlichkeitsfunktion erhielt ich

F(R) = (r^2 (r + 3 R1))/((R1 - R2)^2 (2 R1 + R2))

Nun, und das habe ich dann umgekehrt . Und während draußen Gewitterstürme tobten und der strenge Geruch von Formaldehyd von den eingemachten Organen im Regal aufstieg, erwachte so in meinem Labor beim kalten Lichte der Gasüberschlagslampen das MONSTER . Muahuahuahuahua!!! Sie haben mich wahnsinnig genannt!! Ich werde ihnen zeigen, wer hier wahnsinnig ist!!!!!

Werner Salomon

SeppJ schrieb:

Die Verteilung die ich haben möchte ist
P(R) = R^2/(R2^3/3 - R1^3/3) * Theta(R-R1) * Theta(R2-R)

wo kommt diese Gleichung her. Und was ist die Funktion Theta()? In Statistik bin ich leider ziemlich schwach.

Gruß
Werner