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