001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017 018package org.apache.commons.rng.sampling; 019 020import org.apache.commons.rng.UniformRandomProvider; 021import org.apache.commons.rng.sampling.distribution.NormalizedGaussianSampler; 022import org.apache.commons.rng.sampling.distribution.ZigguratNormalizedGaussianSampler; 023 024/** 025 * Generate vectors <a href="http://mathworld.wolfram.com/SpherePointPicking.html"> 026 * isotropically located on the surface of a sphere</a>. 027 * 028 * <p>Sampling uses:</p> 029 * 030 * <ul> 031 * <li>{@link UniformRandomProvider#nextLong()} 032 * <li>{@link UniformRandomProvider#nextDouble()} 033 * </ul> 034 * 035 * @since 1.1 036 */ 037public class UnitSphereSampler implements SharedStateSampler<UnitSphereSampler> { 038 /** Sampler used for generating the individual components of the vectors. */ 039 private final NormalizedGaussianSampler sampler; 040 /** Space dimension. */ 041 private final int dimension; 042 043 /** 044 * @param dimension Space dimension. 045 * @param rng Generator for the individual components of the vectors. 046 * A shallow copy will be stored in this instance. 047 * @throws IllegalArgumentException If {@code dimension <= 0} 048 */ 049 public UnitSphereSampler(int dimension, 050 UniformRandomProvider rng) { 051 if (dimension <= 0) { 052 throw new IllegalArgumentException("Dimension must be strictly positive"); 053 } 054 055 this.dimension = dimension; 056 sampler = new ZigguratNormalizedGaussianSampler(rng); 057 } 058 059 /** 060 * @param rng Generator for the individual components of the vectors. 061 * @param source Source to copy. 062 */ 063 private UnitSphereSampler(UniformRandomProvider rng, 064 UnitSphereSampler source) { 065 // The Gaussian sampler has no shared state so create a new instance 066 sampler = new ZigguratNormalizedGaussianSampler(rng); 067 dimension = source.dimension; 068 } 069 070 /** 071 * @return a random normalized Cartesian vector. 072 */ 073 public double[] nextVector() { 074 final double[] v = new double[dimension]; 075 076 // Pick a point by choosing a standard Gaussian for each element, 077 // and then normalize to unit length. 078 double normSq = 0; 079 for (int i = 0; i < dimension; i++) { 080 final double comp = sampler.sample(); 081 v[i] = comp; 082 normSq += comp * comp; 083 } 084 085 if (normSq == 0) { 086 // Zero-norm vector is discarded. 087 // Using recursion as it is highly unlikely to generate more 088 // than a few such vectors. It also protects against infinite 089 // loop (in case a buggy generator is used), by eventually 090 // raising a "StackOverflowError". 091 return nextVector(); 092 } 093 094 final double f = 1 / Math.sqrt(normSq); 095 for (int i = 0; i < dimension; i++) { 096 v[i] *= f; 097 } 098 099 return v; 100 } 101 102 /** 103 * {@inheritDoc} 104 * 105 * @since 1.3 106 */ 107 @Override 108 public UnitSphereSampler withUniformRandomProvider(UniformRandomProvider rng) { 109 return new UnitSphereSampler(rng, this); 110 } 111}