// author Stuart Errol Anderson;
// 2015
// stuart.errol.anderson@gmail.com
// n-gon equilateral polygons , angles multiples pi/n
#include <iostream>
#include <stdlib.h> 
#include <fstream>
#include <vector>
#include <cmath>
#include <string>
#include <sstream>
#include <initializer_list>
#define PI 3.14159265358979323846

/////////////////////////////////
using namespace std;
struct Point
{
    double x;
    double y;
};
 
struct Edge
{
	Point p;
	Point q;
};
// Given three colinear points p, q, r, the function checks if
// point q lies on line segment 'pr'
bool onSegment(Point p, Point q, Point r)
{
    if (q.x <= std::max(p.x, r.x) && q.x >= std::min(p.x, r.x) &&
        q.y <= std::max(p.y, r.y) && q.y >= std::min(p.y, r.y))
       return true;
 
    return false;
}

// To find orientation of ordered triplet (p, q, r).
// The function returns following values
// 0 --> p, q and r are colinear
// 1 --> Clockwise
// 2 --> Counterclockwise
int orientation(Point p, Point q, Point r)
{
    // See 10th slides from following link for derivation of the formula
    // http://www.dcs.gla.ac.uk/~pat/52233/slides/Geometry1x1.pdf
    double val = (q.y - p.y) * (r.x - q.x) - (q.x - p.x) * (r.y - q.y);
 	double epsilon = 0.000001;
    if (abs(val) < epsilon) return 0;  // colinear
    return (val > 0)? 1: 2; // clock or counterclock wise
}
 
// The main function that returns true if line segment 'p1q1'
// and 'p2q2' intersect.
bool doIntersect(Point p1, Point q1, Point p2, Point q2)
{
    // Find the four orientations needed for general and
    // special cases
    int o1 = orientation(p1, q1, p2);
    int o2 = orientation(p1, q1, q2);
    int o3 = orientation(p2, q2, p1);
    int o4 = orientation(p2, q2, q1);

    // General case
    if (o1 != o2 && o3 != o4)
        return true;

    // Special Cases
   // p1, q1 and p2 are colinear and p2 lies on segment p1q1
    if (o1 == 0 && onSegment(p1, p2, q1)) return true;

    // p1, q1 and p2 are colinear and q2 lies on segment p1q1
    if (o2 == 0 && onSegment(p1, q2, q1)) return true;

    // p2, q2 and p1 are colinear and p1 lies on segment p2q2
    if (o3 == 0 && onSegment(p2, p1, q2)) return true;

     // p2, q2 and q1 are colinear and q1 lies on segment p2q2
    if (o4 == 0 && onSegment(p2, q1, q2)) return true;

    return false; // Doesn't fall in any of the above cases
}
 
void polygonSearch(std::vector<int>& ixs, const std::vector<int>& endPerIndex, unsigned int currentIndex, unsigned int n) {
    if (currentIndex == ixs.size()) {
        // calculate the x,y components of the vertices of the polygon path; polygon is equilateral; add vectors head to tail
	    unsigned int sum = 0;
		double sx=0.0;
		double sy=0.0;
		vector<double> Vx;
		vector<double> Vy;
		double d =0.0;
		double epsilon = 0.000001;
		string str;
		string filename;
		stringstream ss;
		ofstream myfile;

		for (unsigned int c = 0 ; c != n-1 ; c++) { // first vertex is always (1,0)
			Vx.push_back(cos((c*n-sum)*PI/n));
			sx += Vx[c];
			Vx[c] = sx;
			Vy.push_back(sin((c*n-sum)*PI/n));
			sy += Vy[c];
			Vy[c] = sy;
			sum += ixs[c];   // ixs[c] contains the value of the c-th index.
		}	
		
		// the nth loop is not needed as we can use the fact the sum of the interior angles of an n-gon = (n-2)*180, 		
		unsigned int li = (n-2)*n-sum;  // last index li = (n-2)*n - sum of first n-1 angles [n angle multiples = 180 = pi]	
		if ((sum >= (n-2)*n)||(li > 2*n-1)) return; // reject these cases (get -ve on last index) or angles >= 2pi

		Vx.push_back(cos((n*(n-1)-sum)*PI/n));
		sx += Vx[n-1];
		Vx[n-1] = sx;
		Vy.push_back(sin((n*(n-1)-sum)*PI/n));
		sy += Vy[n-1];
		Vy[n-1] = sy;
		sum += li;   // ixs[c] contains the value of the c-th index.
		d = sqrt(sx*sx+sy*sy);
		if(d<epsilon)  // if the path ends at the origin it is polygon candidate, epsilon can be adjusted to allow a varying error
		{	
			//calculate a canonical ordering by cycling the angle array, converting the array to a number, then picking the cycle order that produces the smallest number
			vector<unsigned int> angle ;
			for (unsigned int j = 0;j!= n-1 ; j++) { 
				angle.push_back(ixs[j]);
			}
			angle.push_back(li);
			int temp = 0;
			long number ;
			unsigned int index1=0;
			int min = pow(n,n);
			for (unsigned int k = 0;k!=n; k++) { 
				//right shift cycle array elements 1 position
				temp = angle[0];
				for (unsigned int j = 0;j!=n; j++) { 
					angle[j] = angle[j+1];
				}
				angle[angle.size()-1] = temp;
				//make a base 2n-2 number from the array values
				number = 0;
				for (unsigned int j = 0;j!=n; j++) { 
					number += angle[j]*pow(n,n-j-1);
				}
				//if it's the smallest so far, make minimum, save shift value
				if (number < min ) { 
					min = number;
					index1 = k;
				}	
				number = 0;
			}
			//right cycle the array values by b places to get canonical angle ordering
			for (unsigned int k = 0;k!=index1+1; k++) { 
				//right shift cycle array elements 1 position
				temp = angle[0];
				for (unsigned int j = 0;j!=n; j++) { 
					angle[j] = angle[j+1];
				}
				angle[angle.size()-1] = temp;
			}
			for (unsigned int k = 0;k!=n; k++) { 
				cout<<angle[k]<<" ";
			}
			cout<<'\n';
			//intersection testing
			//polygon edges meet head to tail, we dont want to test 2 adjacent edges meeting at a vertex, (see special cases in doIntersect definition)
			//intersections can be treated as over/under edges, creating a knot, we count the number of intersections in each solution for classification.
			//paths with no intersections are polygons (knot = false)
			unsigned int intersect_count = 0;
			bool intersect = false;
  			struct Point p1 ; struct Point q1 ;
			struct Point p2 ; struct Point q2 ;
			struct Edge e ;
			vector<Edge> edges;
			for (unsigned int a=0; a<n-1; a++) {
				p1 = {Vx[a], Vy[a]};
			   	q1 = {Vx[a+1],Vy[a+1]};
				e  = {p1, q1}; 
				edges.push_back(e);
			}
			p1 = {Vx[n-1], Vy[n-1]};
		   	q1 = {Vx[0],Vy[0]};
			e   = {p1, q1}; 
			edges.push_back(e);			
			for (unsigned int a=0; a<n; a++) {
				for (unsigned int b = 0; b<n; b++) {
				// detect intersection of edge a with edge b
					if ((a!=b)&&((a+1)!=b)&&((b+1)!=a)&&!((a==(n-1))&&(b==0))&&!((b==(n-1))&&(a==0))) {
						p1 = edges[a].p; q1 = edges[a].q; p2 = edges[b].p; q2 = edges[b].q;						 
						if (doIntersect(p1, q1, p2, q2) == true) {
							intersect_count++;
						}
					}
				}
			}
			//detect concave and convex non-intersecting
			//a convex polygon by definition has no interior angle greater than Pi radians, a polygon which is not convex is concave;
			//a concave polygon has at least one interior angle greater than n = 180 degrees / Pi radians.
			bool concave = false;
			for (unsigned int a=0; a<n; a++) {
				if (angle[a]>n)  concave = true  ;
			}
			if (intersect_count>0)  {
				intersect = true;
			}
			//center the polygon in the A4 page, 595 x 842, center; (297.5,421)
			// a polygon has a maximum diameter n/2 ; cannot extend n/2 in any direction  (flatten a polygon)

			double xscale = 330/(n/2);
			//create the postscript file for the polygon solution - label; knot, concave or convex
			str="%!PS-Adobe- \n";
			str+="% equilateral polygons with pi/"; ss <<n<<" angle multiples\n"; str+=ss.str(); ss.str("");
			str+="% Stuart Anderson Aug 2015\n";
			ss << 297.5 <<" ";
			str+=ss.str(); ss.str("");
			ss << 421 << " translate\n";
			str+=ss.str();ss.str("");
			str+=".02 setlinewidth\n";
			ss << xscale <<" "<< xscale <<" scale\n";
			str+=ss.str();
			ss.str("");
			str+="-0.5 -0.5 moveto\n";
			str+="   /length 1 def\n";
			str+="   /rep "; ss << n <<" def\n"; str+=ss.str(); ss.str("");
			str+="   /rot 180 rep div def\n";
			for (unsigned int c = 0 ; c != n-1 ; c++) {
				str+="	length 0 rlineto\n"; // draw polygon edge 
				ss <<"	"<<n<<" "<<angle[c]<<" sub rot mul rotate\n";  // calculate exterior angle from interior angle Pi * ixs[c] / n and rotate exterior angle
				str+=ss.str();
				ss.str("");
			}
			str+="	length 0 rlineto\n"; // draw polygon edge
			ss<<"% close angle; "<<angle[n-1]<<" sum;"<<sum<<'\n';
			str+=ss.str();
			ss.str("");
			str+="closepath\n";
			str+="gsave\n";
			str+="fill\n";
			str+="grestore\n";
			str+="1 0 0 setrgbcolor\n";
			str+="stroke\n";
			str+="showpage\n";
			

			//knot concave, convex
			if (intersect ==  true ) {
				ss <<"itn_"<<n<<"-";
			} else if (concave == true){
				ss <<"cve_"<<n<<"-"; 
			} else {  
				ss <<"cvx_"<<n<<"-";
			}

			filename += ss.str();
			ss.str("");
			for (unsigned int c = 0 ; c != n ; c++) {
				ss << angle[c]<<"_";
			}		
			ss<<".ps";
			filename += ss.str();			
			myfile.open (filename.c_str());
			myfile << str;
			myfile.close();	
			str ="";
			ss.str("");	
			filename ="";	
									
		}
		sx =0.0;
		sy =0.0;
		d =0;	
		} else {
		    for (ixs[currentIndex] = 1 ; ixs[currentIndex] != endPerIndex[currentIndex] ; ixs[currentIndex]++) { // start at 1, dont want 0 interior angles
		        // Recurse for the next level
		        polygonSearch(ixs, endPerIndex, currentIndex+1, n);
		    }
		}
	}


int main( int argc, const char* argv[] )
{	//argv for n
 	unsigned int n = atoi (argv[1]);
	string str;
	string filename;
	stringstream ss;
	ofstream myfile;
	//The setup for calling polygonSearch from the top level requires two vectors - one for the indexes, ixs
    // and one for the number of iterations at each level, endPerIndex. The example below sets up n-1 nested loops, iterating 2n-1 times at each level:

	vector<int> ixs(n-1, 0);
	vector<int> endPerIndex;
	endPerIndex.push_back(n-1);
	for (unsigned int c = 1 ; c != n ; c++) {  
		endPerIndex.push_back(2*n);
	}
	polygonSearch(ixs, endPerIndex, 0, n);
	return(0);
}