00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 #include "paraboloid.h"
00025 #include "paramset.h"
00026 #include "dynload.h"
00027
00028 using namespace lux;
00029
00030
00031 Paraboloid::Paraboloid(const Transform &o2w, bool ro,
00032 float rad, float z0, float z1,
00033 float tm)
00034 : Shape(o2w, ro) {
00035 radius = rad;
00036 zmin = min(z0,z1);
00037 zmax = max(z0,z1);
00038 phiMax = Radians( Clamp( tm, 0.0f, 360.0f ) );
00039 }
00040 BBox Paraboloid::ObjectBound() const {
00041 Point p1 = Point( -radius, -radius, zmin );
00042 Point p2 = Point( radius, radius, zmax );
00043 return BBox( p1, p2 );
00044 }
00045 bool Paraboloid::Intersect(const Ray &r, float *tHit,
00046 DifferentialGeometry *dg) const {
00047 float phi;
00048 Point phit;
00049
00050 Ray ray;
00051 WorldToObject(r, &ray);
00052
00053 float k = zmax/(radius*radius);
00054 float A = k*(ray.d.x * ray.d.x + ray.d.y * ray.d.y);
00055 float B = 2*k*(ray.d.x * ray.o.x + ray.d.y * ray.o.y) -
00056 ray.d.z;
00057 float C = k*(ray.o.x * ray.o.x + ray.o.y * ray.o.y) -
00058 ray.o.z;
00059
00060 float t0, t1;
00061 if (!Quadratic(A, B, C, &t0, &t1))
00062 return false;
00063
00064 if (t0 > ray.maxt || t1 < ray.mint)
00065 return false;
00066 float thit = t0;
00067 if (t0 < ray.mint) {
00068 thit = t1;
00069 if (thit > ray.maxt) return false;
00070 }
00071
00072 phit = ray(thit);
00073 phi = atan2f(phit.y, phit.x);
00074 if (phi < 0.) phi += 2.f*M_PI;
00075
00076 if (phit.z < zmin || phit.z > zmax || phi > phiMax) {
00077 if (thit == t1) return false;
00078 thit = t1;
00079 if (t1 > ray.maxt) return false;
00080
00081 phit = ray(thit);
00082 phi = atan2f(phit.y, phit.x);
00083 if (phi < 0.) phi += 2.f*M_PI;
00084 if (phit.z < zmin || phit.z > zmax || phi > phiMax)
00085 return false;
00086 }
00087
00088 float u = phi / phiMax;
00089 float v = (phit.z-zmin) / (zmax-zmin);
00090
00091 Vector dpdu(-phiMax * phit.y, phiMax * phit.x, 0.);
00092 Vector dpdv = (zmax - zmin) *
00093 Vector(phit.x / (2.f * phit.z), phit.y / (2.f * phit.z), 1.);
00094
00095 Vector d2Pduu = -phiMax * phiMax *
00096 Vector(phit.x, phit.y, 0);
00097 Vector d2Pduv = (zmax - zmin) * phiMax *
00098 Vector(-phit.y / (2.f * phit.z),
00099 phit.x / (2.f * phit.z),
00100 0);
00101 Vector d2Pdvv = -(zmax - zmin) * (zmax - zmin) *
00102 Vector(phit.x/(4.f*phit.z*phit.z),
00103 phit.y/(4.f*phit.z*phit.z),
00104 0.);
00105
00106 float E = Dot(dpdu, dpdu);
00107 float F = Dot(dpdu, dpdv);
00108 float G = Dot(dpdv, dpdv);
00109 Vector N = Normalize(Cross(dpdu, dpdv));
00110 float e = Dot(N, d2Pduu);
00111 float f = Dot(N, d2Pduv);
00112 float g = Dot(N, d2Pdvv);
00113
00114 float invEGF2 = 1.f / (E*G - F*F);
00115 Normal dndu((f*F - e*G) * invEGF2 * dpdu +
00116 (e*F - f*E) * invEGF2 * dpdv);
00117 Normal dndv((g*F - f*G) * invEGF2 * dpdu +
00118 (f*F - g*E) * invEGF2 * dpdv);
00119
00120 *dg = DifferentialGeometry(ObjectToWorld(phit),
00121 ObjectToWorld(dpdu),
00122 ObjectToWorld(dpdv),
00123 ObjectToWorld(dndu),
00124 ObjectToWorld(dndv),
00125 u, v, this);
00126
00127 *tHit = thit;
00128 return true;
00129 }
00130
00131 bool Paraboloid::IntersectP(const Ray &r) const {
00132 float phi;
00133 Point phit;
00134
00135 Ray ray;
00136 WorldToObject(r, &ray);
00137
00138 float k = zmax/(radius*radius);
00139 float A = k*(ray.d.x * ray.d.x + ray.d.y * ray.d.y);
00140 float B = 2*k*(ray.d.x * ray.o.x + ray.d.y * ray.o.y) -
00141 ray.d.z;
00142 float C = k*(ray.o.x * ray.o.x + ray.o.y * ray.o.y) -
00143 ray.o.z;
00144
00145 float t0, t1;
00146 if (!Quadratic(A, B, C, &t0, &t1))
00147 return false;
00148
00149 if (t0 > ray.maxt || t1 < ray.mint)
00150 return false;
00151 float thit = t0;
00152 if (t0 < ray.mint) {
00153 thit = t1;
00154 if (thit > ray.maxt) return false;
00155 }
00156
00157 phit = ray(thit);
00158 phi = atan2f(phit.y, phit.x);
00159 if (phi < 0.) phi += 2.f*M_PI;
00160
00161 if (phit.z < zmin || phit.z > zmax || phi > phiMax) {
00162 if (thit == t1) return false;
00163 thit = t1;
00164 if (t1 > ray.maxt) return false;
00165
00166 phit = ray(thit);
00167 phi = atan2f(phit.y, phit.x);
00168 if (phi < 0.) phi += 2.f*M_PI;
00169 if (phit.z < zmin || phit.z > zmax || phi > phiMax)
00170 return false;
00171 }
00172 return true;
00173 }
00174 float Paraboloid::Area() const {
00175 return phiMax/12.0f *
00176 (powf(1+4*zmin, 1.5f) - powf(1+4*zmax, 1.5f));
00177 }
00178 Shape* Paraboloid::CreateShape(const Transform &o2w,
00179 bool reverseOrientation, const ParamSet ¶ms) {
00180 float radius = params.FindOneFloat( "radius", 1 );
00181 float zmin = params.FindOneFloat( "zmin", 0 );
00182 float zmax = params.FindOneFloat( "zmax", 1 );
00183 float phimax = params.FindOneFloat( "phimax", 360 );
00184 return new Paraboloid(o2w, reverseOrientation, radius, zmin, zmax, phimax);
00185 }
00186
00187 static DynamicLoader::RegisterShape<Paraboloid> r("paraboloid");