#include "ParticleStream.h"
#include "vObject.h"

#pragma region static

Texture
	*D3D11ParticleStream::def_tex				= NULL;

ID3D11InputLayout *D3D11ParticleStream::IL_Particle		= NULL;

ID3D11VertexShader
	*D3D11ParticleStream::VS_Particle			= NULL;

ID3D11GeometryShader
	*D3D11ParticleStream::GS_Particle			= NULL;

ID3D11PixelShader
	*D3D11ParticleStream::PS_Particle_D			= NULL,
	*D3D11ParticleStream::PS_Particle_E			= NULL,
	*D3D11ParticleStream::PS_Particle_S			= NULL;

ID3D11Buffer
	*D3D11ParticleStream::VB					= NULL,
	*D3D11ParticleStream::cb_GS_Particle_0		= NULL;

GSCB_Particle D3D11ParticleStream::GSCB_Particle_0;

bool D3D11ParticleStream::bShadows = true;

void D3D11ParticleStream::GlobalInit() {
	HRESULT hr;
	ID3DBlob *SBlob = NULL, *EBlob = NULL;

	D3D11_INPUT_ELEMENT_DESC ldesc[ ] = {
		{ "POSITION", 0, DXGI_FORMAT_R32G32B32_FLOAT, 0, 0, D3D11_INPUT_PER_VERTEX_DATA, 0 },
		{ "TEXCOORD", 0, DXGI_FORMAT_R32_FLOAT, 0, 12, D3D11_INPUT_PER_VERTEX_DATA, 0 },
		{ "TEXCOORD", 1, DXGI_FORMAT_R32_FLOAT, 0, 16, D3D11_INPUT_PER_VERTEX_DATA, 0 },
		{ "TEXCOORD", 2, DXGI_FORMAT_R32_UINT, 0, 20, D3D11_INPUT_PER_VERTEX_DATA, 0 },
	};

	CompileFromFile( "Modules\\D3D11Shaders\\ParticleStream.fx", NULL, NULL, "VS_Particle", vs_ver, ShaderCompileFlag, 0, NULL, &SBlob, &EBlob, &hr );
	ShowShaderCompilationError( hr, EBlob );
	HR( Dev->CreateInputLayout( ldesc, 4, SBlob->GetBufferPointer(), SBlob->GetBufferSize(), &IL_Particle ) );
	HR( Dev->CreateVertexShader( SBlob->GetBufferPointer(), SBlob->GetBufferSize(), NULL, &VS_Particle ) );
	REL( SBlob );
	REL( EBlob );

	CompileFromFile( "Modules\\D3D11Shaders\\ParticleStream.fx", NULL, NULL, "GS_Particle", gs_ver, ShaderCompileFlag, 0, NULL, &SBlob, &EBlob, &hr );
	ShowShaderCompilationError( hr, EBlob );
	HR( Dev->CreateGeometryShader( SBlob->GetBufferPointer(), SBlob->GetBufferSize(), NULL, &GS_Particle ) );
	REL( SBlob );
	REL( EBlob );

	CompileFromFile( "Modules\\D3D11Shaders\\ParticleStream.fx", NULL, NULL, "PS_Particle_D", ps_ver, ShaderCompileFlag, 0, NULL, &SBlob, &EBlob, &hr );
	ShowShaderCompilationError( hr, EBlob );
	HR( Dev->CreatePixelShader( SBlob->GetBufferPointer(), SBlob->GetBufferSize(), NULL, &PS_Particle_D ) );
	REL( SBlob );
	REL( EBlob );

	CompileFromFile( "Modules\\D3D11Shaders\\ParticleStream.fx", NULL, NULL, "PS_Particle_E", ps_ver, ShaderCompileFlag, 0, NULL, &SBlob, &EBlob, &hr );
	ShowShaderCompilationError( hr, EBlob );
	HR( Dev->CreatePixelShader( SBlob->GetBufferPointer(), SBlob->GetBufferSize(), NULL, &PS_Particle_E ) );
	REL( SBlob );
	REL( EBlob );

	CompileFromFile( "Modules\\D3D11Shaders\\ParticleStream.fx", NULL, NULL, "PS_Particle_S", ps_ver, ShaderCompileFlag, 0, NULL, &SBlob, &EBlob, &hr );
	ShowShaderCompilationError( hr, EBlob );
	HR( Dev->CreatePixelShader( SBlob->GetBufferPointer(), SBlob->GetBufferSize(), NULL, &PS_Particle_S ) );
	REL( SBlob );
	REL( EBlob );
	D3D11_BUFFER_DESC bdesc;
	ZeroMemory( &bdesc, sizeof(bdesc) );
	bdesc.BindFlags = D3D11_BIND_VERTEX_BUFFER;
	bdesc.ByteWidth = sizeof(PARTICLE_VERTEX)*MAXPARTICLE;
	bdesc.CPUAccessFlags = D3D11_CPU_ACCESS_WRITE;
	bdesc.MiscFlags = 0;
	bdesc.StructureByteStride = 0;
	bdesc.Usage = D3D11_USAGE_DYNAMIC;

	HR( Dev->CreateBuffer( &bdesc, NULL, &VB ) );

	ZeroMemory( &bdesc, sizeof(bdesc) );
	bdesc.BindFlags = D3D11_BIND_CONSTANT_BUFFER;
	bdesc.ByteWidth = 96;
	bdesc.CPUAccessFlags = 0;
	bdesc.MiscFlags = 0;
	bdesc.StructureByteStride = 0;
	bdesc.Usage = D3D11_USAGE_DEFAULT;

	HR( Dev->CreateBuffer( &bdesc, NULL, &cb_GS_Particle_0 ) );

	def_tex = TM->LoadTextureFromFile( "contrail1.dds", 0 );
}

void D3D11ParticleStream::GlobalExit() {
	REL( VS_Particle );
	REL( GS_Particle );
	REL( PS_Particle_D );
	REL( PS_Particle_E );
	REL( PS_Particle_S );
	REL( VB );
	REL( cb_GS_Particle_0 );
}
#pragma endregion

static PARTICLESTREAMSPEC def_pss  = {
	0,								//flags
	8.0,							//creation size
	0.5,							//creation rate
	100.0,							//emission velocity
	0.3,							//velocity randomisation
	8.0,							//lifetime
	0.5,							//growth rate
	3.0,							//atmospheric slowdown
	PARTICLESTREAMSPEC::DIFFUSE,	// render lighting method
	PARTICLESTREAMSPEC::LVL_SQRT,	// mapping from level to alpha
	0.0,							// lmin and lmax levels for mapping
	1.0,
	PARTICLESTREAMSPEC::ATM_PLOG,	// mapping from atmosphere to alpha
	1e-4,							// amin and amax densities for mapping
	1.0,
	NULL
};

D3D11ParticleStream::D3D11ParticleStream( PARTICLESTREAMSPEC *_pss, OBJHANDLE _obj ) : oapi::ParticleStream( gc, _pss ) {
	PARTICLESTREAMSPEC *pss = _pss ? _pss : &def_pss;

	PS.interval = 0.1;
	PS.exp_rate = RAND_MAX/pss->lifetime;
	PS.ipht2 = 0.5/pss->lifetime;

	PS.size0 = pss->srcsize;
	PS.speed = pss->v0;
	PS.vrand = pss->srcspread;
	PS.growthrate = pss->srcspread;
	PS.atmslowdown = pss->atmslowdown;
	PS.pdensity = pss->srcrate;
	diffuse = pss->ltype == PARTICLESTREAMSPEC::DIFFUSE;
	PS.lmap = pss->levelmap;
	PS.lmin = pss->lmin;
	PS.lmax = pss->lmax;
	PS.amap = pss->atmsmap;
	PS.amin = pss->amin;

	switch( PS.amap ) {
	case PARTICLESTREAMSPEC::ATM_PLIN:
		PS.afac = 1.0/(pss->amax - PS.amin);
		break;
	case PARTICLESTREAMSPEC::ATM_PLOG:
		PS.afac = 1.0/log(pss->amax/PS.amin);
		break;
	}

	PS.tex = (TM->Find( (Texture*)pss->tex ) || TM->FindGlobal( (Texture*)pss->tex )) && pss->tex ? (Texture*)pss->tex : def_tex;

	PS.t0 = oapiGetSimTime();
	pfirst = plast = NULL;
	np = 0;

	hVsl = _obj;
	VSL = oapiGetVesselInterface( hVsl );
	vobj = SC->GetVessel( hVsl );
	hSrf = NULL;
	vSrf = NULL;
}

D3D11ParticleStream::~D3D11ParticleStream() {
	while( pfirst ) {
		ParticleSpec *tmp = pfirst;
		pfirst = pfirst->next;
		delete tmp;
	}
}

float D3D11ParticleStream::LevelToAlpha( double level ) {
	switch( PS.lmap ) {
	case PARTICLESTREAMSPEC::LVL_FLAT:
		return PS.lmin;
	case PARTICLESTREAMSPEC::LVL_LIN:
		return (float)level;
	case PARTICLESTREAMSPEC::LVL_SQRT:
		return (float)sqrt(level);
	case PARTICLESTREAMSPEC::LVL_PLIN:
		return max( 0.0f, min( 1.0f, ((float)level - PS.lmin)/(PS.lmax - PS.lmin ) ));
	case PARTICLESTREAMSPEC::LVL_PSQRT:
		return (level <= PS.lmin ? 0.0f : (level >= PS.lmax ? 1.0f : sqrt( (level - PS.lmin)/(PS.lmax - PS.lmin) )));
	}
	return 0.0f;
}

float D3D11ParticleStream::AtmToAlpha( double adens ) {
	switch( PS.amap ) {
	case PARTICLESTREAMSPEC::ATM_FLAT:
		return PS.amin;
	case PARTICLESTREAMSPEC::ATM_PLIN:
		return max( 0.0f, min( 1.0f, (adens - PS.amin)*PS.afac));
	case PARTICLESTREAMSPEC::ATM_PLOG:
		return max( 0.0f, min( 1.0f, log(adens/PS.amin)*PS.afac ));
	}
	return 0.0f;
}

ParticleSpec *D3D11ParticleStream::CreateParticle( const VECTOR3 &pos, const VECTOR3 &vel, double &size, double &alpha ) {
	ParticleSpec *P = new ParticleSpec;

	//set initial params
	P->pos = pos;
	P->vel = vel;
	P->size = (float)size;
	P->alpha0 = (float)alpha;
	P->t0 = (float)simt;
	P->texidx = (rand() & 7);
	P->flag = 0;

	//connection between particles
	P->next = NULL;
	P->prev = plast;
	if( plast )		plast->next = P;
	else			pfirst = P;
	plast = P;
	np++;

	//if there're too many particle - start deleting at the opposite end of the stream
	if( np > MAXPARTICLE )		DeleteParticle( pfirst );

	return P;
}

void D3D11ParticleStream::DeleteParticle( ParticleSpec *P ) {
	//delete from stream
	if( P->prev )	P->prev->next = P->next;
	else			pfirst = P->next;

	if( P->next )	P->next->prev = P->prev;
	else			plast = P->prev;

	//delete
	delete P;
	np--;
}

void D3D11ParticleStream::ComputeSpriteNormal( const VECTOR3 &ppos, D3DXVECTOR3 *vtx ) {
	double uy, uz, vx, vy, vz, len;
	VECTOR3 cdir = unit( ppos );
	if( cdir.y || cdir.z ) {
		uy = cdir.z;
		uz = -cdir.y;
		len = 3.0 / sqrt( uy*uy + uz*uz );
		uy *= len;
		uz *= len;
		vx = cdir.y*cdir.y + cdir.z*cdir.z;
		vy = -cdir.x*cdir.y;
		vz = -cdir.x*cdir.z;
		len = 3.0 / sqrt( vx*vx + vy*vy + vz*vz );
		vx *= len;
		vy *= len;
		vz *= len;
	}
	else {
		uy = 1.0;
		vz = 1.0;
		uz = vx = vy = 0.0;
	}

	static float scale = 0.2294f;
	vtx[0].x = scale*(float)(-cdir.x - vx);
	vtx[0].y = scale*(float)(-cdir.y - uy - vy);
	vtx[0].z = scale*(float)(-cdir.z - uz - vz);

	vtx[1].x = scale*(float)(-cdir.x + vx);
	vtx[1].y = scale*(float)(-cdir.y - uy + vy);
	vtx[1].z = scale*(float)(-cdir.z - uz + vz);

	vtx[2].x = scale*(float)(-cdir.x + vx);
	vtx[2].y = scale*(float)(-cdir.y + uy + vy);
	vtx[2].z = scale*(float)(-cdir.z + uz + vz);

	vtx[3].x = scale*(float)(-cdir.x - vx);
	vtx[3].y = scale*(float)(-cdir.y + uy - vy);
	vtx[3].z = scale*(float)(-cdir.z + uz - vz);
}

void D3D11ParticleStream::Update() {
	ParticleSpec *P, *tmp;

	if( !vobj )
		vobj = SC->GetVessel( hVsl );

	simt = oapiGetSimTime();
	dt = oapiGetSimStep();
	if( np ) {
		if( VSL )	vSrf = SC->GetPlanet( (hSrf = VSL->GetSurfaceRef()) );
		else {
			hSrf = NULL;
			vSrf = NULL;
		}
	}

	//find particles that have to disappear and delete them
	for( P = pfirst; P; ) {
		if( dt*PS.exp_rate > rand() ) {
			tmp = P;
			P = P->next;
			DeleteParticle( tmp );
		}
		else {
			P->pos += P->vel*dt;
			P = P->next;
		}
	}
}

void D3D11ParticleStream::Render() {
	if( !pfirst )	return;	//nothing to render
	if( diffuse )	RenderDiffuse();
	else			RenderEmissive();
}

void D3D11ParticleStream::InitRender() {
	iCtx->IASetPrimitiveTopology( D3D11_PRIMITIVE_TOPOLOGY_POINTLIST );
	iCtx->IASetInputLayout( IL_Particle );
	iCtx->VSSetShader( VS_Particle, NULL, NULL );
	iCtx->GSSetShader( GS_Particle, NULL, NULL );
	iCtx->GSSetConstantBuffers( 0, 1, &cb_GS_Particle_0 );
	iCtx->PSSetConstantBuffers( 0, 1, &cb_D3DXVECTOR4 );
	iCtx->OMSetBlendState( BS_SrcAlpha, D3DXCOLOR( 1, 1, 1, 0 ), 0xFFFFFFFF );
	iCtx->OMSetDepthStencilState( DSS_TestDepthOnly_NoStencil, 0 );
	iCtx->RSSetState( RS_CullBack_Solid );
}

void D3D11ParticleStream::RenderDiffuse() {
	D3DXVECTOR3 SunDir;
	ParticleSpec *P;

	VECTOR3 tpos, cpos = SC->GetCamPos();
	V3toD3( &SunDir, &vobj->GPos );
	D3DXVec3Normalize( &SunDir, &SunDir );

	D3DXVECTOR3 VTX[4];	
	GSCB_Particle_0.VP = *SC->GetVP();
	ComputeSpriteNormal( plast->pos - cpos, VTX );
	GSCB_Particle_0.Light.x = min( 1.0, max( 0.0, -D3DXVec3Dot( &VTX[0], &SunDir ))) + 0.2f;
	GSCB_Particle_0.Light.y = min( 1.0, max( 0.0, -D3DXVec3Dot( &VTX[1], &SunDir ))) + 0.2f;
	GSCB_Particle_0.Light.z = min( 1.0, max( 0.0, -D3DXVec3Dot( &VTX[2], &SunDir ))) + 0.2f;
	GSCB_Particle_0.Light.w = min( 1.0, max( 0.0, -D3DXVec3Dot( &VTX[3], &SunDir ))) + 0.2f;
	GSCB_Particle_0.isShadow = 0;
	iCtx->UpdateSubresource( cb_GS_Particle_0, 0, NULL, &GSCB_Particle_0, 0, 0 );

	iCtx->PSSetShader( PS_Particle_D, NULL, NULL );	
	iCtx->PSSetShaderResources( 0, 1, PS.tex->GetSRV() );
	iCtx->PSSetSamplers( 0, 1, &SS_Linear_Wrap );

	D3D11_MAPPED_SUBRESOURCE SRes;
	memset( &SRes, 0, sizeof(SRes) );
	iCtx->Map( VB, 0, D3D11_MAP_WRITE_DISCARD, 0, &SRes );
	UINT j;
	PARTICLE_VERTEX Vtx, *vtx;
	vtx = (PARTICLE_VERTEX*)SRes.pData;
	for( P = pfirst; P; P = P->next ) {
		tpos = P->pos - cpos;
		Vtx.pos.x = (float)tpos.x, Vtx.pos.y = (float)tpos.y, Vtx.pos.z = (float)tpos.z;
		Vtx.size = (float)P->size;
		Vtx.mix = (float)max( 0.1, P->alpha0*(1.0 - (simt - P->t0)*PS.ipht2 ));
		Vtx.tidx = P->texidx;
		memcpy( vtx, &Vtx, sizeof(PARTICLE_VERTEX) );
		vtx++;
	}
	iCtx->Unmap( VB, 0 );

	const UINT ParticleVertexStride = 24, VBOffset = 0;
	iCtx->IASetVertexBuffers( 0, 1, &VB, &ParticleVertexStride, &VBOffset );

	iCtx->Draw( np, 0 );
}

void D3D11ParticleStream::RenderEmissive() {
	ParticleSpec *P;
	VECTOR3 tpos, cpos = SC->GetCamPos();

	GSCB_Particle_0.VP = *SC->GetVP();
	GSCB_Particle_0.isShadow = 0;
	iCtx->UpdateSubresource( cb_GS_Particle_0, 0, NULL, &GSCB_Particle_0, 0, 0 );

	iCtx->PSSetShader( PS_Particle_E, NULL, NULL );	
	iCtx->PSSetShaderResources( 0, 1, PS.tex->GetSRV() );

	D3D11_MAPPED_SUBRESOURCE SRes;
	memset( &SRes, 0, sizeof(SRes) );
	iCtx->Map( VB, 0, D3D11_MAP_WRITE_DISCARD, 0, &SRes );
	UINT j;
	PARTICLE_VERTEX Vtx, *vtx;
	vtx = (PARTICLE_VERTEX*)SRes.pData;
	for( P = pfirst; P; P = P->next ) {
		tpos = P->pos - cpos;
		Vtx.pos.x = (float)tpos.x, Vtx.pos.y = (float)tpos.y, Vtx.pos.z = (float)tpos.z;
		Vtx.size = (float)P->size;
		Vtx.mix = (float)max( 0.1, P->alpha0*(1.0 - (simt - P->t0)*PS.ipht2 ));
		Vtx.tidx = P->texidx;
		memcpy( vtx, &Vtx, sizeof(PARTICLE_VERTEX) );
		vtx++;
	}
	iCtx->Unmap( VB, 0 );

	if( reentry )
		iCtx->UpdateSubresource( cb_D3DXVECTOR4, 0, NULL, D3DXVECTOR4( 1.0f, 0.7f, 0.5f, 1.0f ), 0, 0 );
	else
		iCtx->UpdateSubresource( cb_D3DXVECTOR4, 0, NULL, D3DXVECTOR4( 1.0f, 1.0f, 1.0f, 1.0f ), 0, 0 );

	const UINT ParticleVertexStride = 24, VBOffset = 0;
	iCtx->IASetVertexBuffers( 0, 1, &VB, &ParticleVertexStride, &VBOffset );

	iCtx->Draw( np, 0 );
}

void D3D11ParticleStream::InitRenderShadow() {
	iCtx->IASetPrimitiveTopology( D3D11_PRIMITIVE_TOPOLOGY_POINTLIST );
	iCtx->IASetInputLayout( IL_Particle );

	iCtx->VSSetShader( VS_Particle, NULL, NULL );
	iCtx->GSSetShader( GS_Particle, NULL, NULL );
	iCtx->PSSetShader( PS_Particle_S, NULL, NULL );

	iCtx->GSSetConstantBuffers( 0, 1, &cb_GS_Particle_0 );
	iCtx->PSSetConstantBuffers( 0, 1, &cb_D3DXVECTOR4 );
	iCtx->OMSetBlendState( BS_SrcAlpha, D3DXCOLOR( 1, 1, 1, 0 ), 0xFFFFFFFF );
	iCtx->OMSetDepthStencilState( DSS_NoDepth_NoStencil, 0 );
	iCtx->RSSetState( RS_CullBack_Solid );
}

void D3D11ParticleStream::RenderShadow() {
	if( !diffuse || !vSrf || !pfirst )	return;

	VECTOR3 tpos, cpos = SC->GetCamPos();
	ParticleSpec *P = pfirst;

	VECTOR3 sd = unit( P->pos );		//shadow projection direction
	VECTOR3 pv0 = P->pos - vSrf->GPos;	//rel. particle pos
	// calculate the intersection of the vessel's shadow with the planet surface
	double fac1 = dotp( sd, pv0 );
	if( fac1 > 0.0 )	return;			// shadow doesn't intersect planet surface
	double arg = fac1*fac1 - (dotp( pv0, pv0 ) - vSrf->rad*vSrf->rad);
	if( arg <= 0.0 )	return;			// shadow doesn't intersect with planet surface
	double a = -fac1 - sqrt( arg );
	VECTOR3 shp = sd*a;					// projection point in global frame
	VECTOR3 hn = unit( shp + pv0 );		// horizon normal in global frame

	GSCB_Particle_0.VP = *SC->GetVP();
	GSCB_Particle_0.isShadow = 1;
	hn *= -1.0;
	GSCB_Particle_0.CDir.x = (float)hn.x;
	GSCB_Particle_0.CDir.y = (float)hn.y;
	GSCB_Particle_0.CDir.z = (float)hn.z;
	iCtx->UpdateSubresource( cb_GS_Particle_0, 0, NULL, &GSCB_Particle_0, 0, 0 );
		
	iCtx->PSSetShaderResources( 0, 1, PS.tex->GetSRV() );

	D3D11_MAPPED_SUBRESOURCE SRes;
	memset( &SRes, 0, sizeof(SRes) );
	iCtx->Map( VB, 0, D3D11_MAP_WRITE_DISCARD, 0, &SRes );
	PARTICLE_VERTEX Vtx, *vtx;
	vtx = (PARTICLE_VERTEX*)SRes.pData;
	DWORD j = 0;
	for( P = pfirst; P; P = P->next ) {
		if( !(P->flag & 1) )
			continue;

		VECTOR3 pvr = P->pos - vSrf->GPos;	// rel. particle position

		double fac1 = dotp( sd, pvr );
		if( fac1 > 0.0 )	break;
		double arg = fac1*fac1 - (dotp( pvr, pvr ) - vSrf->rad*vSrf->rad);
		if( arg <= 0.0 )	break;
		double a = -fac1 - sqrt( arg );

		tpos = P->pos - cpos + sd*a;
		vtx->pos.x = (float)tpos.x, vtx->pos.y = (float)tpos.y, vtx->pos.z = (float)tpos.z;
		vtx->size = (float)P->size;
		vtx->mix = (float)max( 0.1, P->alpha0*(1.0 - (simt - P->t0)*PS.ipht2 ));
		vtx->tidx = P->texidx;
		vtx++;
		j++;
	}
	iCtx->Unmap( VB, 0 );

	const UINT ParticleVertexStride = 24, VBOffset = 0;
	iCtx->IASetVertexBuffers( 0, 1, &VB, &ParticleVertexStride, &VBOffset );

	iCtx->Draw( j, 0 );
}

//==============================================
//				Exhaust stream
//==============================================

ExhaustStream::ExhaustStream( 
		OBJHANDLE _vsobj, const double *lvl, const VECTOR3 *ref,
		const VECTOR3 *dir, PARTICLESTREAMSPEC *pss ) : D3D11ParticleStream( pss, _vsobj )
{
	Attach( hVsl, ref, dir, lvl );
	reentry = false;

}

ExhaustStream::ExhaustStream( 
		OBJHANDLE _vsobj, const double *lvl, const VECTOR3 &ref,
		 const VECTOR3 &dir, PARTICLESTREAMSPEC *pss ) : D3D11ParticleStream( pss, _vsobj )
{
	Attach( hVsl, ref, dir, lvl );
	reentry = false;
}

ExhaustStream::~ExhaustStream() {
}

void ExhaustStream::Update() {
	D3D11ParticleStream::Update();
	
	double alpha0;

	//update particles positions and speed
	if( np ) {
		ParticleSpec *P;
		double
			lng, lat,
			r1,	r2,
			rad;

		UINT j;
		if( vSrf ) {
			//gravitational dv
			VECTOR3 dv = vSrf->GPos - plast->pos;
			double d = length( dv );
			dv *= GGRAV * vSrf->Mass/(d*d*d) * dt;

			//atmospheric dv
			oapiGlobalToEqu( hSrf, pfirst->pos, &lng, &lat, &r1 );
			VECTOR3 av1 = oapiGetWindVector( hSrf, lng, lat, r1 - vSrf->rad, 3 );
			oapiGlobalToEqu( hSrf, plast->pos, &lng, &lat, &r2 );
			VECTOR3 av2 = oapiGetWindVector( hSrf, lng, lat, r2 - vSrf->rad, 3 );
			VECTOR3 dav = (av2 - av1)/np;

			ATMPARAM atmp;
			oapiGetPlanetAtmParams( hSrf, d, &atmp );
			double pref = sqrt(atmp.p)/300.0;
			double slow = exp( -PS.atmslowdown *pref*dt );

			//update speed and position for every particle
			for( P = pfirst, j = 0; P; P = P->next, j++ ) {
				P->vel += dv;
				VECTOR3	av = dav*j + av1;		// atmosphere velocity
				VECTOR3	vv = P->vel - av;		// velocity difference				
				P->vel = vv*slow + av;
				P->size += PS.growthrate*dt;			//growthrate

				//if particle is below surface (reflect a particle)
				VECTOR3 S( P->pos - vSrf->GPos );
				double SLen = length(S);
				if( SLen < vSrf->rad ) {
					VECTOR3 dp = S*(vSrf->rad/SLen - 1.0);
					P->pos += dp;

					static double dv_scale = length(vv)*0.2;
					VECTOR3 dV = {
						((double)rand()/(double)RAND_MAX - 0.5)*dv_scale,
						((double)rand()/(double)RAND_MAX - 0.5)*dv_scale,
						((double)rand()/(double)RAND_MAX - 0.5)*dv_scale	};

					dV += vv;
					normalise(S);
					VECTOR3 vv2 = dV - S*dotp( S, dV );
					double vv2_len = length(vv2); 
					if( vv2_len )
						vv2 *= 0.5*length(vv)/vv2_len;
					vv2 += S*( ((double)rand()/(double)RAND_MAX)*dv_scale);
					P->vel = vv2 + av;
					double rnd = (double)rand()/(double)RAND_MAX;
					P->pos += (vv2 - vv)*dt*rnd;
				}
			}
		}
	}

	//create new particles if needed
	if( level && *level > 0.0 && (alpha0 = LevelToAlpha( *level ) * AtmToAlpha( VSL->GetAtmDensity() )) > 0.0 ) {
		if( simt > PS.t0 + PS.interval ) {
			VECTOR3 vv, vr;
			VSL->GetGlobalVel( vv );
			vr = mul( vobj->mROT, *dir )*(-PS.speed);
			while( simt > PS.t0 + PS.interval ) {
				//create new particle
				double
					dt = simt - PS.t0 - PS.interval,
					dv_scale = PS.speed*PS.vrand;
				VECTOR3 dV = {
					((double)rand()/(double)RAND_MAX - 0.5)*dv_scale,
					((double)rand()/(double)RAND_MAX - 0.5)*dv_scale,
					((double)rand()/(double)RAND_MAX - 0.5)*dv_scale	};

				ParticleSpec *P = CreateParticle( mul( vobj->mROT, *pos ) + vobj->GPos + (vr + dV)*dt, vv + vr + dV, PS.size0, alpha0 );
				P->size += PS.growthrate*dt;

				if( diffuse && vSrf && bShadows ) {
					double lng, lat, alt;
					static const double eps = 1e-2;
					oapiGlobalToEqu( hSrf, P->pos, &lng, &lat, &alt );
					alt -= vSrf->rad;
					if( alt*eps < VSL->GetSize() )
						P->flag |= 1;	//render shadow
				}

				//determine next interval
				PS.t0 += PS.interval;
				if( PS.speed > 10.0 )
					PS.interval = max( 0.015, PS.size0/(PS.pdensity*(0.1*VSL->GetAirspeed() + PS.size0)) );
				else
					PS.interval = 1.0/PS.pdensity;
				PS.interval *= (double)rand()/(double)RAND_MAX + 0.5;
			}
		}
	}
	else
		PS.t0 = simt;
}

//==============================================
//				Reentry stream
//==============================================

ReentryStream::ReentryStream( PARTICLESTREAMSPEC *pss,	OBJHANDLE _vsobj ) : D3D11ParticleStream( pss, _vsobj )
{
	llevel = 1.0;
	Attach( hVsl, _V( 0.0, 0.0, 0.0 ), _V( 0.0, 0.0, 0.0 ), &llevel );
	reentry = true;
}

ReentryStream::~ReentryStream() {
}

void ReentryStream::Update() {
	D3D11ParticleStream::Update();

	if( !hRef )	return;
	double friction = ( VSL ? VSL->GetDynPressure()*VSL->GetAirspeed() : 0.0 );
	double alpha0;
	UINT j;

	if( np ) {
		ParticleSpec *P;
		double lng, lat, r1, r2, rad;
		if( vSrf ) {
			//atmospheric dv
			oapiGlobalToEqu( hSrf, pfirst->pos, &lng, &lat, &r1 );
			VECTOR3 av1 = oapiGetWindVector( hSrf, lng, lat, r1 - vSrf->rad, 3 );
			oapiGlobalToEqu( hSrf, plast->pos, &lng, &lat, &r2 );
			VECTOR3 av2 = oapiGetWindVector( hSrf, lng, lat, r2 - vSrf->rad, 3 );
			VECTOR3 dav = (av2 - av1)/np;

			for( P = pfirst, j = 0; P; P = P->next, j++ ) {
				VECTOR3 av = dav*j + av1;
				VECTOR3 vv = P->vel - av;
				double slow = exp( -PS.atmslowdown*dt );
				P->vel = vv*slow + av;
				P->size += PS.growthrate*dt;
			}
		}
	}

	if( friction > 0.0 && (alpha0 = AtmToAlpha( friction )) > 0.01 ) {
		if( simt > PS.t0 + PS.interval ) {
			VECTOR3 vv, av;
			VSL->GetGlobalVel( vv );

			if( hSrf ) {
				double lng, lat, r;
				oapiGlobalToEqu( hSrf, vobj->GPos, &lng, &lat, &r );
				av = oapiGetWindVector( hSrf, lng, lat, r - vSrf->rad, 3 );
			}
			else
				av = vv;

			while( simt > PS.t0 + PS.interval ) {
				//create new particle
				double _dt = simt - PS.t0 - PS.interval;
				double ebt = exp( -PS.atmslowdown*_dt );
				double dv_scale = VSL->GetAirspeed()*PS.vrand;

				VECTOR3 dv = {
					((double)rand()/(double)RAND_MAX-0.5)*dv_scale,
					((double)rand()/(double)RAND_MAX-0.5)*dv_scale,
					((double)rand()/(double)RAND_MAX-0.5)*dv_scale	};

				VECTOR3 dx = (vv-av)*(1.0-ebt)/PS.atmslowdown + av*_dt;
				CreateParticle( vobj->GPos + dx - vv*_dt, (vv + dv - av)*ebt + av, PS.size0, alpha0 );

				//determine next interval
				PS.t0 += PS.interval;
				PS.interval = max( 0.015, PS.size0/(PS.pdensity*(0.1*VSL->GetAirspeed() + PS.size0)) );
				PS.interval *= (double)rand()/(double)RAND_MAX + 0.5;
			}
		}
	}
	else
		PS.t0 = simt;
}