scene.cpp

蒙特卡罗方法可以有效地解决复杂的工程问题

	m_Primitive[m_Primitives++] = new Primitive( Primitive::TRIANGLE, v[1], v[5], v[2] );
	m_Primitive[m_Primitives++] = new Primitive( Primitive::TRIANGLE, v[5], v[6], v[2] );
	// top plane
	m_Primitive[m_Primitives++] = new Primitive( Primitive::TRIANGLE, v[4], v[5], v[1] );
	m_Primitive[m_Primitives++] = new Primitive( Primitive::TRIANGLE, v[4], v[1], v[0] );
	// bottom plane
	m_Primitive[m_Primitives++] = new Primitive( Primitive::TRIANGLE, v[6], v[7], v[2] );
	m_Primitive[m_Primitives++] = new Primitive( Primitive::TRIANGLE, v[7], v[3], v[2] );
}

void Scene::AddPlane( vector3 a_P1, vector3 a_P2, vector3 a_P3, vector3 a_P4, Material* a_Mat )
{
	Vertex* v[4];
	v[0] = new Vertex( a_P1, 0, 0 );
	v[1] = new Vertex( a_P2, 0, 0 );
	v[2] = new Vertex( a_P3, 0, 0 );
	v[3] = new Vertex( a_P4, 0, 0 );
	m_Primitive[m_Primitives] = new Primitive( Primitive::TRIANGLE, v[0], v[1], v[3] );
	m_Primitive[m_Primitives++]->SetMaterial( a_Mat );
	m_Primitive[m_Primitives] = new Primitive( Primitive::TRIANGLE, v[1], v[2], v[3] );
	m_Primitive[m_Primitives++]->SetMaterial( a_Mat );
}

bool Scene::InitScene( Surface* a_MsgSurf )
{
	Material* mat;
	int x;
	vector3 p1, p2;
	switch (m_State)
	{
	case 0:
		a_MsgSurf->Print( "constructing geometry", 2, 2, 0xffffffff );
		break;
	case 1:
		m_Primitive = new Primitive*[20000];
		m_Primitives = 0;
		m_Light = new Light*[MAXLIGHTS];
		m_Lights = 0;
		// ground plane
		mat = new Material();
		mat->SetParameters( 0.0f, 0.0f, Color( 0.4f, 0.3f, 0.3f ), 1.0f, 0.0f );
		AddPlane( vector3( -13, -4.4f, -5.5f ), vector3( 13, -4.4f, -5.5f ),
				  vector3(  13, -4.4f, 29 ), vector3( -13, -4.4f, 29 ), mat );
		// back plane
		mat = new Material();
		mat->SetParameters( 0.0f, 0.0f, Color( 0.5f, 0.3f, 0.5f ), 0.6f, 0.0f );
		AddPlane( vector3( -13, -4.4f, 8 ), vector3( 13, -4.4f, 16 ),
				  vector3(  13,  7.4f, 16 ), vector3( -13, 7.4f, 8 ), mat );
		// ceiling plane
		mat = new Material();
		mat->SetParameters( 0.0f, 0.0f, Color( 0.4f, 0.7f, 0.7f ), 0.5f, 0.0f );
		AddPlane( vector3(  13, 7.4f, -5.5f ), vector3( -13, 7.4f, -5.5f ),
				  vector3( -13, 7.4f, 29 ), vector3( 13, 7.4f, 29 ), mat );
		// big sphere
		// m_Primitive[m_Primitives] = new Primitive( Primitive::SPHERE, vector3( 2, 0.8f, 3 ), 2.5f );
		// m_Primitive[m_Primitives++]->GetMaterial()->SetParameters( 0.2f, 0.8f, Color( 0.7f, 0.7f, 1.0f ), 0.2f, 0.8f );
		// small sphere
		m_Primitive[m_Primitives] = new Primitive( Primitive::SPHERE, vector3( -5.5f, -0.5, 7 ), 2 );
		m_Primitive[m_Primitives++]->GetMaterial()->SetParameters( 0.5f, 0.0f, Color( 0.7f, 0.7f, 1.0f ), 0.1f, 0.8f );
	#if 0
		// area lights
		m_Light[m_Lights++] = new Light( Light::AREA, vector3( -1, 6, 4 ), vector3( 1, 6, 4 ), vector3( -1, 6, 6 ), Color( 0.7f, 0.7f, 0.7f ) );
		m_Light[m_Lights++] = new Light( Light::AREA, vector3( -1, 6, -1 ), vector3( 1, 6, -1 ), vector3( -1, 6, 1 ), Color( 0.7f, 0.7f, 0.7f ) );
	#else
		// point lights
		m_Light[m_Lights++] = new Light( Light::POINT, vector3( 0, 5, 5 ), Color( 0.4f, 0.4f, 0.4f ) );
		m_Light[m_Lights++] = new Light( Light::POINT, vector3( -3, 5, 1 ), Color( 0.6f, 0.6f, 0.8f ) );
	#endif
		// extra sphere
		m_Primitive[m_Primitives] = new Primitive( Primitive::SPHERE, vector3( -1.5f, -3.8f, 1 ), 1.5f );
		m_Primitive[m_Primitives++]->GetMaterial()->SetParameters( 0.0f, 0.8f, Color( 1.0f, 0.4f, 0.4f ), 0.2f, 0.8f );
		// grid
		for ( x = 0; x < 8; x++ ) for ( int y = 0; y < 7; y++ )
		{
			m_Primitive[m_Primitives] = new Primitive( Primitive::SPHERE, vector3( -4.5f + x * 1.5f, -4.3f + y * 1.5f, 10 ), 0.3f );
			m_Primitive[m_Primitives++]->GetMaterial()->SetParameters( 0.0f, 0.0f, Color( 0.3f, 1.0f, 0.4f ), 0.6f, 0.6f );
		}
		for ( x = 0; x < 8; x++ ) for ( int y = 0; y < 8; y++ )
		{
			m_Primitive[m_Primitives] = new Primitive( Primitive::SPHERE, vector3( -4.5f + x * 1.5f, -4.3f, 10.0f - y * 1.5f ), 0.3f );
			m_Primitive[m_Primitives++]->GetMaterial()->SetParameters( 0.0f, 0.0f, Color( 0.3f, 1.0f, 0.4f ), 0.6f, 0.6f );
		}
		for ( x = 0; x < 16; x++ ) for ( int y = 0; y < 8; y++ )
		{
			m_Primitive[m_Primitives] = new Primitive( Primitive::SPHERE, vector3( -8.5f + x * 1.5f, 4.3f, 10.0f - y ), 0.3f );
			m_Primitive[m_Primitives++]->GetMaterial()->SetParameters( 0.0f, 0.0f, Color( 0.3f, 1.0f, 0.4f ), 0.6f, 0.6f );
		}
		mat = new Material();
		mat->SetParameters( 0.9f, 0, Color( 0.9f, 0.9f, 1 ), 0.3f, 0.7f );
		mat->SetRefrIndex( 1.3f );
		// mat->SetTexture( new Texture( "textures/wood.tga" ) );
		Load3DS( "meshes/knot.3ds", mat, vector3( 0, 0.5f, 4 ), 6 );
		break;
	case 2:
		a_MsgSurf->Print( "building kdtree - please wait", 2, 10, 0xffffffff );
		break;
	case 3:
		// build the kd-tree
		p1 = vector3( -14, -6, -6 ), p2 = vector3( 14, 8, 30 );
		m_Extends = aabb( p1, p2 - p1 );
		m_KdTree = new KdTree();
		m_KdTree->Build( this );
		break;
	default:
		return true;
	};
	m_State++;
	return false;
}

// -----------------------------------------------------------
// KdTree class implementation
// -----------------------------------------------------------
KdTree::KdTree()
{
	m_Root = new KdTreeNode();
}

void KdTree::Build( Scene* a_Scene )
{
	for ( int p = 0; p < a_Scene->GetNrPrimitives(); p++ )
		m_Root->Add( a_Scene->GetPrimitive( p ) );
	int prims = a_Scene->GetNrPrimitives();
	aabb sbox = a_Scene->GetExtends();
	m_SPool = new SplitList[prims * 2 + 8];
	for ( int i = 0; i < (prims * 2 + 6); i++ ) m_SPool[i].next = &m_SPool[i + 1];
	m_SPool[i].next = 0;
	m_SList = 0;
	Subdivide( m_Root, sbox, 0, prims );
}

void KdTree::InsertSplitPos( real a_SplitPos )
{
	// insert a split position candidate in the list if unique
	SplitList* entry = m_SPool;
	m_SPool = m_SPool->next;
	entry->next = 0;
	entry->splitpos = a_SplitPos;
	entry->n1count = 0;
	entry->n2count = 0;
	if (!m_SList) m_SList = entry; else
	{
		if (a_SplitPos < m_SList->splitpos)
		{
			entry->next = m_SList;
			m_SList = entry;
		}
		else if (a_SplitPos == m_SList->splitpos)
		{
			entry->next = m_SPool; // redundant; recycle
			m_SPool = entry;
		}
		else
		{
			SplitList* list = m_SList;
			while ((list->next) && (a_SplitPos >= list->next->splitpos)) 
			{
				if (a_SplitPos == list->next->splitpos)
				{
					entry->next = m_SPool; // redundant; recycle
					m_SPool = entry;
					return;
				}
				list = list->next;
			}
			entry->next = list->next;
			list->next = entry;
		}
	}
}

void KdTree::Subdivide( KdTreeNode* a_Node, aabb& a_Box, int a_Depth, int a_Prims )
{
	// recycle used split list nodes
	if (m_SList)
	{
		SplitList* list = m_SList;
		while (list->next) list = list->next;
		list->next = m_SPool;
		m_SPool = m_SList, m_SList = 0;
	}
	// determine split axis
	vector3 s = a_Box.GetSize();
	if ((s.x >= s.y) && (s.x >= s.z)) a_Node->SetAxis( 0 );
	else if ((s.y >= s.x) && (s.y >= s.z)) a_Node->SetAxis( 1 );
	int axis = a_Node->GetAxis();
	// make a list of the split position candidates
	ObjectList* l = a_Node->GetList();
	real p1, p2;
	real pos1 = a_Box.GetPos().cell[axis];
	real pos2 = a_Box.GetPos().cell[axis] + a_Box.GetSize().cell[axis];
	bool* pright = new bool[a_Prims];
	float* eleft = new float[a_Prims], *eright = new float[a_Prims];
	Primitive** parray = new Primitive*[a_Prims];
	int aidx = 0;
	while (l)
	{
		Primitive* p = parray[aidx] = l->GetPrimitive();
		pright[aidx] = true;
		p->CalculateRange( eleft[aidx], eright[aidx], axis );
		aidx++;
		if (p->GetType() == Primitive::SPHERE)
		{
			p1 = p->GetCentre().cell[axis] - p->GetRadius();
			p2 = p->GetCentre().cell[axis] + p->GetRadius();
			if ((p1 >= pos1) && (p1 <= pos2)) InsertSplitPos( p1 );
			if ((p2 >= pos1) && (p2 <= pos2)) InsertSplitPos( p2 );
		}
		else
		{
			for ( int i = 0; i < 3; i++ )
			{
				p1 = p->GetVertex( i )->GetPos().cell[axis];
				if ((p1 >= pos1) && (p1 <= pos2)) InsertSplitPos( p1 );
			}
		}
		l = l->GetNext();
	}
	// determine n1count / n2count for each split position
	aabb b1, b2, b3 = a_Box, b4 = a_Box;
	SplitList* splist = m_SList;
	float b3p1 = b3.GetPos().cell[axis];
	float b4p2 = b4.GetPos().cell[axis] + b4.GetSize().cell[axis];
	while (splist)
	{
		b4.GetPos().cell[axis] = splist->splitpos;
		b4.GetSize().cell[axis] = pos2 - splist->splitpos;
		b3.GetSize().cell[axis] = splist->splitpos - pos1;
		float b3p2 = b3.GetPos().cell[axis] + b3.GetSize().cell[axis];
		float b4p1 = b4.GetPos().cell[axis];
		for ( int i = 0; i < a_Prims; i++ ) if (pright[i])
		{
			Primitive* p = parray[i];
			if ((eleft[i] <= b3p2) && (eright[i] >= b3p1))
				if (p->IntersectBox( b3 )) splist->n1count++;
			if ((eleft[i] <= b4p2) && (eright[i] >= b4p1))
				if (p->IntersectBox( b4 )) splist->n2count++; else pright[i] = false;
		}
		else splist->n1count++;
		splist = splist->next;
	}
	delete pright;
	// calculate surface area for current node
	real SAV = 0.5f / (a_Box.w() * a_Box.d() + a_Box.w() * a_Box.h() + a_Box.d() * a_Box.h());
	// calculate cost for not splitting
	real Cleaf = a_Prims * 1.0f;
	// determine optimal split plane position
	splist = m_SList;
	real lowcost = 10000;
	real bestpos = 0;
	while (splist)
	{
		// calculate child node extends
		b4.GetPos().cell[axis] = splist->splitpos;
		b4.GetSize().cell[axis] = pos2 - splist->splitpos;
		b3.GetSize().cell[axis] = splist->splitpos - pos1;
		// calculate child node cost
		real SA1 = 2 * (b3.w() * b3.d() + b3.w() * b3.h() + b3.d() * b3.h());
		real SA2 = 2 * (b4.w() * b4.d() + b4.w() * b4.h() + b4.d() * b4.h());
		real splitcost = 0.3f + 1.0f * (SA1 * SAV * splist->n1count + SA2 * SAV * splist->n2count);
		// update best cost tracking variables
		if (splitcost < lowcost)
		{
			lowcost = splitcost;
			bestpos = splist->splitpos;
			b1 = b3, b2 = b4;
		}
		splist = splist->next;
	}
	if (lowcost > Cleaf) return;
	a_Node->SetSplitPos( bestpos );
	// construct child nodes
	KdTreeNode* left = s_MManager->NewKdTreeNodePair();
	int n1count = 0, n2count = 0, total = 0;
	// assign primitives to both sides
	float b1p1 = b1.GetPos().cell[axis];
	float b2p2 = b2.GetPos().cell[axis] + b2.GetSize().cell[axis];
	float b1p2 = b1.GetPos().cell[axis] + b1.GetSize().cell[axis];
	float b2p1 = b2.GetPos().cell[axis];
	for ( int i = 0; i < a_Prims; i++ )
	{
		Primitive* p = parray[i];
		total++;
		if ((eleft[i] <= b1p2) && (eright[i] >= b1p1)) if (p->IntersectBox( b1 )) 
		{
			left->Add( p );
			n1count++;
		}
		if ((eleft[i] <= b2p2) && (eright[i] >= b2p1)) if (p->IntersectBox( b2 )) 
		{
			(left + 1)->Add( p );
			n2count++;
		}
	}
	delete eleft;
	delete eright;
	delete parray;
	s_MManager->FreeObjectList( a_Node->GetList() );
	a_Node->SetLeft( left );
	a_Node->SetLeaf( false );
	if (a_Depth < MAXTREEDEPTH)
	{
		if (n1count > 2) Subdivide( left, b1, a_Depth + 1, n1count );
		if (n2count > 2) Subdivide( left + 1, b2, a_Depth + 1, n2count );
	}
}

// -----------------------------------------------------------
// KdTreeNode class implementation
// -----------------------------------------------------------

void KdTreeNode::Add( Primitive* a_Prim )
{
	ObjectList* lnode = KdTree::s_MManager->NewObjectList();
	lnode->SetPrimitive( a_Prim );
	lnode->SetNext( GetList() );
	SetList( lnode );
}

}; // namespace Raytracer