I'm trying to implement a particle system (in 3D) were an emitter can have a direction and spread.

The direction is defined as a 3D normalized vector.

The spread is a scalar value (0.0 - 1.0) which represents the range, relative to the direction, that particles can start moving.

• So a spread=0.0 would be the same as the direction vector. This special case makes the spread value moot.
• A spread=0.5 would be a half-sphere, with the dome-side oriented in the same direction as the direction vector.
• A spread=1.0 would be a whole sphere. This special case makes the direction vector moot.

The following example illustrates the problem in 2D. The red coloring highlights the spread.

In each picture (except the last one) the direction vector is going straight up (x=0, y=1).

The last picture (parden my poor drawing) is supposed to be the normalized vector (x=0.7, y=0.7).

[attachment=33650:example.png]

#include <iostream>
#include <cmath>
#include <ctime>
 
float const Pi = std::atan(1.0f) * 4.0f;
 
float uniform_random(float a, float b) {
  return a + drand48() * (b - a);
}
 
struct Vector {
  float x, y, z;
  
  Vector(float x, float y, float z) : x(x), y(y), z(z) {
  }
};
 
std::ostream &operator<<(std::ostream &os, Vector const &v) {
  return os << '(' << v.x << ',' << v.y << ',' << v.z << ')';
}
 
Vector operator+(Vector const &v1, Vector const &v2) {
  return Vector(v1.x + v2.x, v1.y + v2.y, v1.z + v2.z);
}
 
Vector operator*(float s, Vector const &v) {
  return Vector(s * v.x, s * v.y, s * v.z);
}
 
float dot(Vector const &v1, Vector const &v2) {
  return v1.x * v2.x + v1.y * v2.y + v1.z * v2.z;
}
 
float length(Vector const &v) {
  return std::sqrt(dot(v, v));
}
 
Vector cross_product(Vector const &v1, Vector const &v2) {
  return Vector(v1.y * v2.z - v1.z * v2.y,
                v1.z * v2.x - v1.x * v2.z,
                v1.x * v2.y - v1.y * v2.x);
}
 
Vector normalize(Vector const &v) {
  return (1.0f / length(v)) * v;
}
 
Vector random_direction(Vector direction, float spread) {
  // Make an orthogonal basis whose third vector is along `direction'
  Vector b3 = normalize(direction);
  Vector different = (std::abs(b3.x) < 0.5f) ? Vector(1.0f, 0.0f, 0.0f) : Vector(0.0f, 1.0f, 0.0f);
  Vector b1 = normalize(cross_product(b3, different));
  Vector b2 = cross_product(b1, b3);
 
  // Pick (x,y,z) randomly around (0,0,1)
  float z = uniform_random(std::cos(spread * Pi), 1);
  float r = std::sqrt(1.0f - z * z);
  float theta = uniform_random(-Pi, +Pi);
  float x = r * std::cos(theta);
  float y = r * std::sin(theta);
 
  // Construct the vector that has coordinates (x,y,z) in the basis formed by b1, b2, b3
  return x * b1 + y * b2 + z * b3;
}
 
int main() {
  srand48(std::time(0));
  
  for (int i = 0; i < 10; ++i)
    std::cout << random_direction(Vector(1.0f, 0.0f, 0.0f), 0.1f) << '\n';
}

Pick a random number z uniformly between cos(spread*pi) and 1. Pick a random angle theta. Make x = sqrt(1-z^2)*cos(theta), y = sqrt(1-z^2)*sin(theta). Now apply to (x,y,z) a rotation that turns (0,0,1) into the desired direction. Voilà!

That makes sense.

Some code for you (I was bored). The meat is in the function `random_direction'.

Thank you so much! You made my day.

You should be aware that this distribution is not uniform and if you emit a lot of particles you are likely to see clustering, especially toward the poles. If the clustering is acceptable then this is probably the best solution, otherwise you may prefer to use a function which generates uniformly random points on the surface of a sphere and iterates until it finds a point that lands inside of your cone. Normalizing the output of a multivariate Gaussian distribution is a good way to get a uniform distribution on the surface of a sphere, see https://en.m.wikipedia.org/wiki/Multivariate_normal_distribution

You must be new here... ;)

You must be new here... ;)

Relatively new, yup... is it that obvious? :)

Alright, I did some investigation and I retract my earlier statement. I was aware that uniformly sampling polar angles would create bunching near the poles and I didn't expect uniformly sampling along a fixed axis to be any better, since it was not immediately obvious to me that slicing a sphere at uniform intervals along the z axis would divide the sphere surface into pieces of equal area. It looks like http://mathworld.wolfram.com/SpherePointPicking.html agrees that this is actually the case. I will edit the above post accordingly.

Edit: looks like your edit beat me to the punch :)