/* (C) 2012 Claudio Rocchini - CC BY 3.0 */  #include <stdio.h> #include <stdlib.h> #include <math.h> #include <vector> #include <algorithm>  const double PI = 3.1415926535897932384626433832795; inline double drand() { return double(rand())/RAND_MAX; }  class point2 { public: 	double x,y; 	point2() {} 	point2( double nx, double ny ): x(nx),y(ny) {} };  class seed { public: 	point2 p;		// seed position 	double a;		// main angle 	double l[2];	// extents };  class inter {		// intersection public: 	int i,j; 	double ti; double tj; 	inline inter(); 	inline inter( int ni, int nj, double nti, double ntj ) : i(ni),j(nj),ti(nti),tj(ntj) {} 	inline bool operator< ( const inter & ii ) const { return fabs(ti)<fabs(ii.ti); } };  bool intersects( const seed & s1, const seed & s2, double & t1, double & t2 ) { 	double dx1 = cos(s1.a); double dy1 = sin(s1.a); 	double dx2 = cos(s2.a); double dy2 = sin(s2.a); 	const double de = dx1*dy2-dy1*dx2; 	if(fabs(de)<1e-8) return false; 	double ox1 = s1.p.x; double oy1 = s1.p.y; 	double ox2 = s2.p.x; double oy2 = s2.p.y; 	t1 = (-ox1*dy2+ox2*dy2+dx2*oy1-dx2*oy2)/de; 	t2 = ( dx1*oy1-dx1*oy2-dy1*ox1+dy1*ox2)/de; 	return t1<s1.l[0] && t1>-s1.l[1] && t2<s2.l[0] && t2>-s2.l[1]; }  int main() { 	const double SX = 800; const double SY = 800; const int N = 256; 	int i,j; 	static seed seeds[N]; 	std::vector<inter> inters; std::vector<inter>::iterator ii; 	 	for(i=0;i<N;++i) { 		seeds[i].p.x = drand()*SX; seeds[i].p.y = drand()*SY; 		seeds[i].a   = PI*drand(); // ortho version: drand()>= 0.5 ? 0 : PI/2;  		seeds[i].l[0] = SX*10; seeds[i].l[1] = SX*10; 	} 	for(i=0;i<N-1;++i) for(j=i+1;j<N;++j) { 		double t1,t2; 		if(intersects(seeds[i],seeds[j],t1,t2)) { 			inters.push_back( inter(i,j,t1,t2) ); 			inters.push_back( inter(j,i,t2,t1) ); 	}	} 	std::sort(inters.begin(),inters.end()); 	for(ii=inters.begin();ii!=inters.end();++ii) {	// compute intersection on time 		if(ii->ti < seeds[ii->i].l[0] && ii->ti > -seeds[ii->i].l[1] && 	       ii->tj < seeds[ii->j].l[0] && ii->tj > -seeds[ii->j].l[1]) { 			if( fabs(ii->ti) > fabs(ii->tj) ) { 				if(ii->ti>0) seeds[ii->i].l[0] =  ii->ti; 				else         seeds[ii->i].l[1] = -ii->ti; 	}   }   } 	 	for(i=0;i<N;++i) {	// page clipping 		seed s; 		s.p.x = 0; s.p.y = 0; s.a = 0; s.l[0] = SX; s.l[1] = 0; 		double t1,t2; 		if(intersects(seeds[i],s,t1,t2)) { if(t1>0) seeds[i].l[0] = t1; else seeds[i].l[1] = -t1; } 		s.p.y = SY; 		if(intersects(seeds[i],s,t1,t2)) { if(t1>0) seeds[i].l[0] = t1; else seeds[i].l[1] = -t1; } 		s.p.x = 0; s.p.y = 0; s.l[0] = SY; s.l[1] = 0; s.a = PI/2; 		if(intersects(seeds[i],s,t1,t2)) { if(t1>0) seeds[i].l[0] = t1; else seeds[i].l[1] = -t1; } 		s.p.x = SX; 		if(intersects(seeds[i],s,t1,t2)) { if(t1>0) seeds[i].l[0] = t1; else seeds[i].l[1] = -t1; } 	} 	FILE * fo = fopen("gilbert.svg","w");	// output the results 	fprintf(fo, 		"<?xml version=\"1.0\" encoding=\"UTF-8\" standalone=\"no\"?>\n" 		"<svg xmlns=\"http://www.w3.org/2000/svg\" version=\"1.0\" width=\"%g\" height=\"%g\" id=\"%s\">\n" 		,SX,SY,"girlbert" 	); 	fprintf(fo,"<g style=\"stroke:#000000;stroke-width:2;stroke-linecap:round;fill:none\">\n"); 	for(i=0;i<N;++i){ 		fprintf(fo,"<line x1=\"%6.2f\" y1=\"%6.2f\" x2=\"%6.2f\" y2=\"%6.2f\"/>\n" 			,seeds[i].p.x+seeds[i].l[0]*cos(seeds[i].a) 			,seeds[i].p.y+seeds[i].l[0]*sin(seeds[i].a) 			,seeds[i].p.x+seeds[i].l[1]*cos(seeds[i].a+PI) 			,seeds[i].p.y+seeds[i].l[1]*sin(seeds[i].a+PI) 		); 	} 	fprintf(fo,"</g>\n"); 	fprintf(fo,"</svg>\n"); 	fclose(fo); 	return 0; }