Convex Hull এবং গ্রাহামের স্ক্যান

কনভেক্স হাল ……

নেটে সার্চ দিলাম । এই বিষয়ের উপর বাংলায়  টিউটোরিয়াল পেলাম । তাই নিজেই লিখতে বসে গেলাম বুঝার জন্য....ভুল হলে ধরিয়ে দিবেন আশা করি...।।

ফার্স্ট অফ অল… হোয়াট ইজ কনভেক্স , অর কনকেভ ?

ইন্টেরস্টেড যারা তারা এই লিঙ্ক এ গিয়ে এর উত্তর ভাল ভাবে পেতে পারেন ।

কনভেক্স না কনকেভ:

ইন এ নাটশেল …

২টা শর্ত পূরণ হলে কোন পলিগন কে কনভেক্স বলা যেতে পারে ।

# যদি প্রত্যেকটা ইন্টার্নাল অ্যাঙ্গেল এর মান ১৮০°বা তার কম হয় ।

# যেকোন ২ টা ভার্টিস বা শীর্ষবিন্দু এর মধ্যে যতগুলো লাইন সেগমেন্ট টানা যাবে সবগুলোই ওই পলিগনের ভিতরে বা বর্ডারে উপর থাকবে ।

এই ২টা শর্ত ঠিকঠাক ভাবে পূরণ হলেই কেবল মাত্র আমরা যেকোন পলিগন কে কনভেক্স বললে পারি, আদারওয়াইজ ওইটা কনকেভ ।

সেকেন্ড অফ অল… হোয়াট ইজ কনভেক্স হাল … (Convex HULL)

যদি তোমার কাছে একগাদা পয়েন্ট থাকে , তাহলে এই সবগুলো পয়েন্ট কে কভার করে মিনিমাম এরিয়ার যে কনভেক্স শেপ টা তুমি কল্পনা করতে পারবে সেইটাই কনভেক্স হাল ।

মানে এমন একটা শেপ তোমাকে কল্পনা করতে হবে যাতে সবগুলো পয়েন্ট ওই শেপ এর ভিতরে, না হয় ওই শেপ এর বাউন্ডারিতে থাকে , বাইরে নয় , অবশ্যই সেই শেপ টা কনভেক্স হতে হবে , তা না হলেতো সেইটা কনভেক্স হাল হবে না :P

তো আমরা জানলাম কনভেক্স কি জিনিস , আর ও জানলাম কনভেক্স হাল কি জিনিস ।

এখন তোমাকে যদি আমি একগাদা পয়েন্ট 2-D ডাইমেনশন এ দিয়ে দেই ,আর যদি বলি তোমাকে সেই পয়েন্ট গুলার কনভেক্স হাল বের করতে হবে অর্থাৎ যত গুলা পয়েন্ট তোমাকে দেয়া হইছে তার মধ্যে কোন কোন পয়েন্ট গুলা কে ভার্টেক্স বা শীর্ষবিন্দু হিসেবে সিলেক্ট করে পলিগন বানালে সেইটা কনভেক্স হাল হবে। তাহলে তুমি কি করবে ?

বাংলা সোজা কথা …

ইনপুট…

এক গাদা পয়েন্ট

আউটপুট

এক গাদা পয়েন্ট

এখন অ্যালগোরিদম …

ফার্স্ট কথা , কনভেক্স হাল বের করার জন্য অনেক অ্যালগোরিদম আছে , এখানে (n log n) রান করে এমন একটা পপুলার অ্যালগোরিদম নিয়েই আলোচনা করব ।

অ্যালগোরিদম টার নাম গ্রাহামের স্ক্যান …

এখন অ্যালগোরিদম টা শিখতে গেলে বা ভালোভাবে বুঝতে গেলে আরো কিছু বেসিক জ্যামিতির জিনিস তোমাকে জানতে হবে ।

১. কাউন্টার ক্লকওয়াইজ টেস্ট(তিনটা পয়েন্ট কাউন্টার ক্লক ওয়াইজ এ আছে নাকি ক্লক ওয়াইজ এ সেটা টেস্টিং)

তুমি যদি রেখার ঢাল সম্পর্কে জানো তাহলে এইটা ঢাল হিসেব করেই টেস্ট করা যায়। যেমন মনে করো A(x1,y1) ,B(x2,y2) ,C(x3,y3) তিনটা পয়েন্ট, তোমাকে টেস্ট করতে হবে ABC কাউন্টার ক্লকওয়াইজ ।

হিসাব করার সুবিধার জন্য প্রথমে ধরে নেই A পয়েন্টা সবচেয়ে বাম দিকে x co-ordinate এর দিক দিয়ে আছে , অর্থাৎ x1 < x2 এবং x1 < x3।

এখন ABC কাউন্টার ক্লক ওয়াইজ হবে যদি AB এর ঢাল AC এর ঢাল থেকে কম হয় ।

এখন AB এর ঢাল (y2-y1)/(x2-x1) , AC এর ঢাল (y3-y1)/(x3-x1)

বা, (y2-y1)/(x2-x1) < (y3-y1)/(x3-x1)

বা, (y2-y1)*(x3-x1) < (x2-x1)*(y3-y1)

যদি এই ইনইকুয়েলিটি হোল্ড করে তাহলে কাউন্টার ক্লক ওয়াইজ আর না হলে ক্লক ওয়াইজ।

আবার এই একই টেস্ট,ভেক্টরের সাহায্যেও করা যায় ।

আমরা ভেক্টরের ক্রস প্রডাক্ট বের করে নির্ধারণ করতে পারি এই টেস্টটা।

AB x AC যদি পজিটিভ রেসাল্ট দেয় তাহলে কাউন্টার ক্লক ওয়াইজ আর নেগেটিভ হলে ক্লক ওয়াইজ, আর 0 হলে লিনিয়ার মানে তিনটা বিন্দু ই এক সরলরেখায় অবস্থিত।

অর্থাৎ

| (x2-x1)……(y2-y1) |

| (x3-x1)……(y3-y1) |

এই ম্যাট্রিক্সের ডিটারমিনেন্ট 0 থেকে বেশি হতে হবে।

এখন আমরা জানি CCW test মানে কাউন্টার ক্লক ওয়াইজ টেস্ট ।

এইবার তোমাকে জানতে হবে একটা বেসিক ডাটা স্ট্রাকচার সম্পর্কে , আর সেটা হলো STACK

Stack হল LIFO , LAST IN FIRST OUT.

এইটা তুমি এইভাবে চিন্তা করতে পার, তোমার বাসায় তোমার মা , বা কাজের লোক যখন প্লেট পরিষ্কার করতে থাকে তখন দেখ একগাদা প্লেট একসাথে একটার উপর একটা সাজিয়ে রাখে , তখন সবার শেষে যে প্লেট টা রাখে সেটা থাকে সবার উপরে , এইটা কে বলে স্ট্যাক (STACK).

এখন তুমি এই স্ট্যাক থেকে কোন কিছু বের করতে চাইলে নিশ্চয়ই সবার উপরের টা আগে বের করবে মানে সবার লাস্টে যেটা স্ট্যাকে পুশ করেছো সেইটা ।

স্ট্যাক সম্পর্কে আরো জানতে চাইলে …

স্ট্যাক

এখন আমরা গ্রাহামের স্ক্যান সম্পর্কে জানব।

১. প্রথমে আমরা bottom most এবং rightmost পয়েন্ট টাকে খুঁজে বের করব। এই পয়েন্টটাকে বলা হয় পিভট। এইবার এই পিভটের সাপেক্ষে সকল পয়েন্টকে সর্ট করতে হবে , উইথ রেস্পেক্ট টু অ্যাঙ্গেল । মানে অন্য পয়েন্ট গুলা পিভট এর সাথে ছোট কোণ করেছে যারা তারা আগে থাকবে, যদি টাই হয়ে যায় (মানে ২ টা পয়েন্ট যদি একই কোণ করে), তাহলে কাছাকাছি যে পয়েন্ট টা থাকবে (পিভটের), সেইটা আগে আসবে।

২.এইবার আমাদের একটা স্ট্যাক মেইন্টেইন করতে হবে , স্ট্যাকে এই পয়েন্ট গুলা সবই ঢুকবে , কিন্তু সব অপারেশনের শেষে যে পয়েন্ট গুলা কনভেক্স হাল গঠন করে না , মানে শীর্ষবিন্দু ছাড়া অন্য বিন্দুগুলা পপ আউট করবে মানে, স্ট্যাক থেকে বের হয়ে যাবে। গ্রেইট… কিন্তু , আমরা বুঝব কিভাবে , কোন শীর্ষবিন্দু কনভেক্স হাল গঠন করবে কি করবে না ? কিছু কি মনে পড়ছে? আমরা টিউটোরিয়ালটার শুরুতে দেখে এসেছি কাউন্টার ক্লক ওয়াইজ টেস্ট সম্পর্কে।

৩.আমাদের স্ট্যাকে সবার শুরুতে থাকবে ২ টা পয়েন্ট,পিভট আর সর্ট করার পর সবার শেষ বিন্দু,এইবার এই ২ টা বিন্দুকে সোয়াপ করব(কেন করব? তা না হলে প্রথমেই কাউন্টার ক্লক ওয়াইজ টেস্ট ফেইল করবে) এরপর সাথে আরেকটা বিন্দু যোগ করব(৩য় বিন্দু টি হচ্ছে সর্ট করার পর প্রথম বিন্দু) যদি পাশাপাশি ৩ টা বিন্দু(অ্যাঙ্গেলের ভিত্তিতে সর্ট করা , মনে রেখো)কাউন্টার ক্লক ওয়াইজ টেস্ট পাস করে তাহলে ওই বিন্দুগুলা কনভেক্স হাল এর অন্তর্ভুক্ত, আর যদি ফেইল করে , তাহলে ২য় বিন্দুটিকে পপ করে দিবো, এরপর নেক্সট পয়েন্ট টা কে প্রসেস করার জন্য স্টাকে ঢুকাব।

এক্সপ্লানেশন

প্রথমে স্ট্যাকএর ভিতরে থাকবে পিভট,আর সর্ট করার পর সবার শেষের পয়েন্ট টা। কাউন্টারক্লক ওয়াইজ টেস্ট যাতে ফেইল না করে এরজন্য এই ২টা পয়েন্টকে সোয়াপ করে নিতে হয় (GrahamScan ফাংশানের while লুপের আগের রুটিন দ্রষ্টব্য) , তাহলে প্রথমে ছিল 0-11, সোয়াপ এর পর হল 11-0, এরপর নেক্সট পয়েন্ট হল 1, তাহলে স্ট্যাকের ভিতরে ২ টা পয়েন্ট আর 1 এই তিনটা পয়েন্ট নিয়ে CCW Test করা হল। রেজাল্ট এল পজিটিভ,তাই 1 স্ট্যাকে পুশ করব। তাহলে স্ট্যাকের চেহারা হবে 11-0-1, এবার নেক্সট পয়েন্ট হল 2 , তাহলে এইবার CCW Test হবে 0-1-2 কে নিয়ে ,এবং এইটাও পজিটিভ রেজাল্ট দিবে, তাই 2 কেও পুশ করব।

এখন স্ট্যাক টা দাঁড়াবে 11-0-1-2,নেক্সট পয়েন্ট 3 ,তাই এইবার CCW test হবে 1-2-3 কে নিয়ে । চিত্রে লক্ষ্য করো,এইটা লেফট টার্ন না, রাইট টার্ন, মানে ক্লক ওয়াইজ। সো 2 কে পপ করে দিবো। সো স্ট্যাকের চেহারা হবে 11-0-1-3….

এইভাবে চলতে থাকবে ,সব পয়েন্ট ঘুরে লাস্টে এসে আবার 11 তে এসে থামবে।

অবশেষে যা থাকবে সেই পয়েন্ট গুলাই হবে কনভেক্স হালের শীর্ষবিন্দু,যদিও শেষের পয়েন্ট টা ডুপ্লিকেট হবে, কারণ পলিগন একটা বিন্দু থেকে শুরু হয়ে আবার ওই বিন্দুতেই শেষ হয়।

এই হচ্ছে গ্রাহামের স্ক্যানের মূল বিষয়গুলো…এখন কোড 

typedef struct POINT{

    int x,y;

}point;

 point pivot;

vector<point>polygon,HULL;    

int turn(point p, point q, point r) {

    int result = (r.x - q.x) * (p.y - q.y) - (r.y - q.y) * (p.x - q.x);

    if (result < 0) return -1;  // P,Q,R is a right turn

    if (result > 0) return 1;  // P,Q,R is a left turn

    return 0; // P,Q,R is a straight line, i.e. P, Q, R are collinear

}

bool ccw(point p, point q, point r) {

    return (turn(p, q, r) > 0); 

 }

int area2(point a,point b,point c){

    return a.x * b.y - a.y * b.x + b.x * c.y - b.y * c.x + c.x * a.y - c.y * a.x;

}//assuming all points are integer type

int dist2(point a,point b){

    int dx = a.x-b.x;

    int dy = a.y-b.y;

    return dx*dx+dy*dy;     //note no need to compute sqrt, why?

}

bool angle_cmp(point a,point b){

    if(area2(pivot,a,b)==0){    //collinear

        return dist2(pivot,a) < dist2(pivot,b)   //which one is closer

    }

     

    int d1x = a.x - pivot.x, d1y = a.y - pivot.y;

    int d2x = b.x - pivot.x, d2y = b.y - pivot.y;

     

    return (atan2((double)d1y,(double)d1x) - atan2((double)d2y,(double)d2x)) < 0;

}

vector <point> GrahamScan(vector<point>Polygon){   

//first of all we need to find P0, pivot with lowest Y, if tie with rightmost X

     

    int i,P0=0,N = Polygon.size();

    for(i=1;i<N;i++){

        if((Polygon[i].y < Polygon[P0].y)||

        (Polygon[i].y==Polygon[P0].y  && Polygon[i].x > Polygon[P0].x)){

            P0 = i;

        }

    }

    //swap selected vertex with Polygon[0]

    point temp = Polygon[0];

    Polygon[0] = Polygon[P0];

    Polygon[P0] = temp;

    //now Polygon[0] contains the pivot

    //sort w.r.t.angle skipping Polygon[0]

    pivot = Polygon[0];

    sort(++Polygon.begin(),Polygon.end(),angle_cmp);

    //ccw test

    stack <point> S;

    point prev,now;

    S.push(Polygon[N-1]); //put 2 starting vertices in the stack

    S.push(Polygon[0]);

    //start checking the rest

    i = 1;

    while(i<N){

        now = S.top();

        S.pop();

        prev = S.top();

        S.push(now);    //trick to get 2nd item from top of S

     

    

        if(ccw(prev,now,Polygon[i])){

            S.push(Polugon[i]);  //if these 3 points make a left turn

            i++;                 //accept

        }else {

     

            S.pop(); //pop this point until we have a left turn

        }

     

    }

    vector <point>HULL;

    while(!S.empty()){    //from stack

        HULL.push_back(S.top());

        S.pop();

    }

     

    HULL.pop_back()//the last one is the duplicate of the first element

    return HULL;

     

}

Convex Fence Code:

#include<stdio.h>

#include<math.h>

#include<algorithm>

#define EPS 1e-10

// +0.0000000001

#define S(x) ((x)*(x))

#define Z(x) (fabs(x) < EPS)

using namespace std;

struct Vector{

double x,y;

    void scan(){scanf("%lf %lf",&x,&y);}

void print(){printf("(%lf,%lf)\n",x,y);}

Vector(double x=0,double y=0){

this->x=x;

this->y=y;

}

double mag2(){

return S(x)+S(y);

}

double mag(){return sqrt(mag2());}

};

double crossp(Vector a, Vector b){

return a.x*b.y - a.y*b.x;

}

Vector operator-(Vector a,Vector b){

return Vector(a.x-b.x, a.y-b.y); }

int Turn(Vector &V0,Vector &V1,Vector &V2)

{

    double v=crossp((V1-V0),(V2-V0));

    if(fabs(v)<EPS)return 0;

    if(v>0)return 1;

    return -1;

}

Vector V[50007];

int hull[50007];

bool cmp(Vector V1,Vector V2)

{

    int t=Turn(V[0],V1,V2);

    if(!t)

    {

        return (V1-V[0]).mag2()<(V2-V[0]).mag2();

    }

    return t>0;

}

int cHull(int N)

{

    int i,id=0;

    int top=0;

    if(N<=2)

    {

        hull[top++]=0;

        if(N==2)hull[top++]=1;

        return top;

    }

    for(i=1;i<N;i++)

    {

        if(V[i].y==V[id].y)

        {

            if(V[i].x<V[id].x)id=i;

        }

        if(V[i].y<V[id].y)id=i;

    }

    swap(V[0],V[id]);

    sort(V+1,V+N,cmp);

    hull[top++]=0;

    hull[top++]=1;

    for(i=2;i<N;i++)

    {

        while(Turn(V[hull[top-2]],V[hull[top-1]],V[i])<1&&top>1)top--;

        hull[top++]=i;

    }

    return top;

}

int main()

{

   // freopen("in.txt","r",stdin);

    int T,N,i,ks=0;

    double r;

    scanf("%d",&T);

    while(T--)

    {

        ks++;

        scanf("%d%lf",&N,&r);

        for(i=0;i<N;i++)V[i].scan();

        int top=cHull(N);

        hull[top]=hull[0];

        double ans=0;

        for(i=0;i<top;i++)

        {

            ans+=(V[hull[i]]-V[hull[i+1]]).mag();

            V[hull[i]].print();

        }

        ans+=2*acos(-1)*r;

        printf("Case %d: %.10lf\n",ks,ans);

    }

    return 0;

}

End program.................By IM@N