//Jacek Matulewski, e-mail: jacek@fizyka.umk.pl

#ifndef UKLADYPUNKTOWMATERIALNYCH3_H
#define UKLADYPUNKTOWMATERIALNYCH3_H

#include "PunktMaterialny.h"
#include "UkladyPunktowMaterialnych2.h"

#define _USE_MATH_DEFINES
#include <math.h>
#include <time.h>

#include <vector>

class Galaktyka : public ZbiorPunktowMaterialnychZObszaremZabronionym
{
	private:
		Wektor polozenieMin,polozenieMax;
		double G;		
		double progOdleglosci;
		bool podzialNaKomorki;
		bool periodyczneWarunkiBrzegowe;		
		
		Wektor rozmiarKomorki;		
		int ileKomorekX,ileKomorekY,ileKomorekZ,ileKomorek;
		std::vector<int>* komorki;
		int* indeksyKomorek;

	public:
		Galaktyka(int iloscGwiazd,double stalaGrawitacji,Wektor polozenieMin,Wektor polozenieMax,double progOdleglosci)
			:ZbiorPunktowMaterialnychZObszaremZabronionym(iloscGwiazd),
			 G(stalaGrawitacji),progOdleglosci(progOdleglosci),
			 periodyczneWarunkiBrzegowe(true),podzialNaKomorki(true),
			 polozenieMin(polozenieMin),polozenieMax(polozenieMax)
		{
			//obszarZabroniony=new Pudlo(1,0,polozenieMin.X,polozenieMax.X,polozenieMin.Y,polozenieMax.Y,polozenieMin.Z,polozenieMax.Z);
			
			srand((unsigned)time(NULL));
			//UstawGwiazdyLosowo(polozenieMin,polozenieMax,0,1,1,1);
			UstawGwiazdyWObwarzanek(Wektor(-1.5,0,0),1,0.5,Wektor(0,0,1),Wektor(1,0,0));

			/*
			PunktMaterialny* punkt0=PobierzPunktMaterialny(0);
			punkt0->UstawMase(100);
			punkt0->UstawPolozenie(Wektor(0,0,0));
			punkt0->UstawPredkosc(Wektor(0,0,0));
			punkt0->UstawKolor(1,1,0);

			PunktMaterialny* punkt1=PobierzPunktMaterialny(1);
			punkt1->UstawMase(10);
			punkt1->UstawPolozenie(Wektor(-1,0,0));
			punkt1->UstawPredkosc(Wektor(0,0.5,0));
			punkt1->UstawKolor(1,0,0);

			PunktMaterialny* punkt2=PobierzPunktMaterialny(2);
			punkt2->UstawMase(10);
			punkt2->UstawPolozenie(Wektor(1,0,0));
			punkt2->UstawPredkosc(Wektor(0,-1,0));
			punkt2->UstawKolor(0,1,0);
			*/



			if(podzialNaKomorki)
			{
				//generowanie komorek
				ileKomorekX=(int)ceil((polozenieMax.X-polozenieMin.X)/progOdleglosci);
				ileKomorekY=(int)ceil((polozenieMax.Y-polozenieMin.Y)/progOdleglosci);
				ileKomorekZ=(int)ceil((polozenieMax.Z-polozenieMin.Z)/progOdleglosci);
				ileKomorek=ileKomorekX*ileKomorekY*ileKomorekZ;
				rozmiarKomorki.X=(polozenieMax.X-polozenieMin.X)/ileKomorekX;
				rozmiarKomorki.Y=(polozenieMax.Y-polozenieMin.Y)/ileKomorekY;
				rozmiarKomorki.Z=(polozenieMax.Z-polozenieMin.Z)/ileKomorekZ;

				komorki=new std::vector<int>[ileKomorek];
				for(int i=0;i<ileKomorek;++i) komorki[i]=std::vector<int>();
				indeksyKomorek=new int[ilosc];
				PrzydzielPunktyDoKomorek();
			}
		}

		~Galaktyka()
		{
			if(podzialNaKomorki)
			{
				CzyscKomorki();
				delete [] komorki;
				delete [] indeksyKomorek;
			}
		}

	private:
		void CzyscKomorki()
		{
			if(!podzialNaKomorki) return;
			for(int i=0;i<ileKomorek;++i) komorki[i].clear();
		}

		void PrzydzielPunktyDoKomorek()
		{	
			if(!podzialNaKomorki) return;
			CzyscKomorki();
			for(int i=0;i<ilosc;++i)
			{
				PunktMaterialny* punktI=PobierzPunktMaterialny(i);
				Wektor polozenie=punktI->Polozenie();
				int ktoraKomorkaX=(int)floor((polozenie.X-polozenieMin.X)/rozmiarKomorki.X);
				int ktoraKomorkaY=(int)floor((polozenie.Y-polozenieMin.Y)/rozmiarKomorki.Y);
				int ktoraKomorkaZ=(int)floor((polozenie.Z-polozenieMin.Z)/rozmiarKomorki.Z);
				indeksyKomorek[i]=ktoraKomorkaX+ileKomorekX*ktoraKomorkaY+ileKomorekX*ileKomorekY*ktoraKomorkaZ;
				if(indeksyKomorek[i]>=0 && indeksyKomorek[i]<ileKomorek) komorki[indeksyKomorek[i]].push_back(i);
			}
		}

		void UstawGwiazdyLosowo(Wektor polozenieMin,Wektor polozenieMax,double predkoscMin,double predkoscMax,double masaMin,double masaMax)
		{
			for(int i=0;i<ilosc;++i)
			{
				bool wynikTestuNakrywania=false;
				do
				{
					double polozenieX=(double)rand()/(RAND_MAX+1)*(polozenieMax.X-polozenieMin.X)+polozenieMin.X;
					double polozenieY=(double)rand()/(RAND_MAX+1)*(polozenieMax.Y-polozenieMin.Y)+polozenieMin.Y;
					double polozenieZ=(double)rand()/(RAND_MAX+1)*(polozenieMax.Z-polozenieMin.Z)+polozenieMin.Z;				
					PunktMaterialny* punkt=PobierzPunktMaterialny(i);				
					punkt->UstawPolozenie(Wektor(polozenieX,polozenieY,polozenieZ));
					double predkoscKatPhi=(double)rand()/(RAND_MAX+1)*2*M_PI;
					double szybkosc=(double)rand()/(RAND_MAX+1)*(predkoscMax-predkoscMin)+predkoscMin;
					punkt->UstawPredkosc(szybkosc*Wektor(cos(predkoscKatPhi),sin(predkoscKatPhi),0));
					double masa=(double)rand()/(RAND_MAX+1)*(masaMax-masaMin)+masaMin;
					punkt->UstawMase(masa);				
					punkt->UstawKolor(i/(float)ilosc,1-i/(float)ilosc,1);

					//test nakladania punktow
					wynikTestuNakrywania=false;
					for(int j=0;j<i;++j)
					{
						PunktMaterialny* punktJ=PobierzPunktMaterialny(j);
						if((punkt->Polozenie()-punktJ->Polozenie()).Dlugosc()==0) //oslabic warunek?
							wynikTestuNakrywania=true;
					}
				}
				while(wynikTestuNakrywania);
			}						
		}

		void UstawGwiazdyWObwarzanek(Wektor srodekMasy,double promienMin,double grubosc,Wektor normalna,Wektor predkosc)
		{
			for(int i=0;i<ilosc;++i)
			{
				bool wynikTestuNakrywania=false;
				do
				{					
					double odlegloscSloncaOdSrodkaMasy=(double)rand()/(RAND_MAX+1)*(grubosc-promienMin)+promienMin;
					double katPhi=(double)rand()/(RAND_MAX+1)*2*M_PI;
					const double predkoscLiniowa=1.6;
					PunktMaterialny* punkt=PobierzPunktMaterialny(i);
					punkt->UstawPolozenie(srodekMasy+odlegloscSloncaOdSrodkaMasy*Wektor(cos(katPhi),sin(katPhi),0));
					Wektor kierunekPredkosci=(punkt->Polozenie()-srodekMasy)^normalna;
					punkt->UstawPredkosc(predkosc+kierunekPredkosci*predkoscLiniowa);
					punkt->UstawKolor(i/(float)ilosc,1-i/(float)ilosc,1);

					//test nakladania punktow
					wynikTestuNakrywania=false;
					for(int j=0;j<i;++j)
					{
						PunktMaterialny* punktJ=PobierzPunktMaterialny(j);
						if((punkt->Polozenie()-punktJ->Polozenie()).Dlugosc()==0) //oslabic warunek?
							wynikTestuNakrywania=true;
					}
				}
				while(wynikTestuNakrywania);
			}
		}

		Wektor Roznica(Wektor p2,Wektor p1) const
		{
			Wektor roznica=p2-p1;
			if(periodyczneWarunkiBrzegowe)
			{
				roznica.X=(fabs(roznica.X)<0.5*rozmiarKomorki.X)?roznica.X:(roznica.X>0?(roznica.X-rozmiarKomorki.X):(roznica.X+rozmiarKomorki.X));
				roznica.Y=(fabs(roznica.Y)<0.5*rozmiarKomorki.Y)?roznica.Y:(roznica.Y>0?(roznica.Y-rozmiarKomorki.Y):(roznica.Y+rozmiarKomorki.Y));
				roznica.Z=(fabs(roznica.Z)<0.5*rozmiarKomorki.Z)?roznica.Z:(roznica.Z>0?(roznica.Z-rozmiarKomorki.Z):(roznica.Z+rozmiarKomorki.Z));			
			}
			return roznica;
		}

		Wektor Sila(int indeks) const
		{
			//return Wektor(0,0,0);

			PunktMaterialny* punkt=PobierzPunktMaterialny(indeks);
			Wektor sila;

			if(!podzialNaKomorki)
			{								
				for(int i=0;i<ilosc;++i)
					if(i!=indeks) 
					{
						PunktMaterialny* punktI=PobierzPunktMaterialny(i);
						Wektor r12=Roznica(punktI->Polozenie(),punkt->Polozenie());
						//Wektor r12=punktI->Polozenie()-punkt->Polozenie();
						double odleglosc2=SQR(r12.Dlugosc());
						if(odleglosc2<progOdleglosci)
						{
							r12.Normuj();
							sila+=(punkt->Masa()*punktI->Masa()/odleglosc2)*r12;
						}
					}				
			}
			else
			{
				int indeksKomorkiCentralnej=indeksyKomorek[indeks]; //komorka, w ktorej jest gwiazda wskazana przez indeks
				int ix0=indeksKomorkiCentralnej%ileKomorekX;
				int iy0=((indeksKomorkiCentralnej-ix0)/ileKomorekX)%ileKomorekY;
				int iz0=(indeksKomorkiCentralnej-ix0*iy0)/(ileKomorekX*ileKomorekY);

				const int ileSasiadow=0;
				for(int ix=max(0,ix0-ileSasiadow);ix<=min(ix0+ileSasiadow,ileKomorekX-1);++ix)
					for(int iy=max(0,iy0-ileSasiadow);iy<=min(iy0+ileSasiadow,ileKomorekY-1);++iy)
						for(int iz=max(0,iz0-ileSasiadow);iz<=min(iz0+ileSasiadow,ileKomorekZ-1);++iz)
						{
							int indeksKomorki=ix+ileKomorekX*iy+ileKomorekX*ileKomorekY*iz;							
							if(!komorki[indeksKomorki].empty())
							{
								for(int i=0;i<(int)komorki[indeksKomorki].size();++i)
								{
									PunktMaterialny* punktI=PobierzPunktMaterialny(komorki[indeksKomorki][i]);
									Wektor r12=Roznica(punktI->Polozenie(),punkt->Polozenie());
									//Wektor r12=punktI->Polozenie()-punkt->Polozenie();
									double odleglosc2=SQR(r12.Dlugosc());
									if(odleglosc2==0) continue;
									if(odleglosc2<progOdleglosci)
									{
										r12.Normuj();
										sila+=(punkt->Masa()*punktI->Masa()/odleglosc2)*r12;
									}
								}
							}
						}
			}

			return G*sila;
		}		

		void PoPrzygotowaniuRuchu(double krokCzasowy)
		{
			if(periodyczneWarunkiBrzegowe)
			{			
				for(int i=0;i<ilosc;++i)
				{
					PunktMaterialny* punktI=PobierzPunktMaterialny(i);
					Wektor polozenie=punktI->Polozenie();
					bool odswiez=false;
					if(polozenie.X<polozenieMin.X) {polozenie.X+=polozenieMax.X-polozenieMin.X; odswiez=true;}
					if(polozenie.X>polozenieMax.X) {polozenie.X-=polozenieMax.X-polozenieMin.X; odswiez=true;}
					if(polozenie.Y<polozenieMin.Y) {polozenie.Y+=polozenieMax.Y-polozenieMin.Y; odswiez=true;}
					if(polozenie.Y>polozenieMax.Y) {polozenie.Y-=polozenieMax.Y-polozenieMin.Y; odswiez=true;}
					if(polozenie.Z<polozenieMin.Z) {polozenie.Z+=polozenieMax.Z-polozenieMin.Z; odswiez=true;}
					if(polozenie.Z>polozenieMax.Z) {polozenie.Z-=polozenieMax.Z-polozenieMin.Z; odswiez=true;}
					if(odswiez)
					{
						punktI->UstawPolozenie(polozenie);
						punktI->PrzygotujRuch(Sila(i),krokCzasowy,algorytmEulera);
					}
				}
			}
			
			if(podzialNaKomorki) PrzydzielPunktyDoKomorek();
		}
};

#endif