%\smallskip
%\noindent \textbf{GSM and Geometry Code}
\subsection{GSM and Geometry Code}


\lstset{language=C++}
%012345678901234567890123456789012345678901234
\begin{lstlisting}
#include <complex>
#include <queue>
#include <set>
#include <limits>
#include "prelude.cpp"

using namespace std;

typedef complex<double> P;
namespace std {
  // a < b if: a has smaller y value, or,
  // if they have equal y value, smaller x
  bool operator<(const P& a, const P& b) {
    return imag(b-a) > 0. ||
           imag(b-a) == 0. &&
           real(b-a) > 0.; }
}

double e = 1e-9; // roundoff error fuzz

double cross(P a, P b)
  { return -imag(a*conj(b)); }
double   dot(P a, P b)
  { return  real(a*conj(b)); }
P       proj(P v, P u)
  { return u * P(dot(u,v)/dot(u,u)); }

#define angle arg
#define len abs
P unit(P v) { return v/len(v); }

typedef pair<P,P> Seg;

Seg from_pts(P a, P b) {
    return Seg(b-a, a); }

P tail(Seg s) { return s.second; }
P  pos(Seg s) { return s.second; }
P  dir(Seg s) { return s.first; }
P head(Seg s) { return pos(s) + dir(s); }

// angle between two lines
double angle(Seg a, Seg b) {
  return abs(angle(dir(b))
           - angle(dir(a))); }

namespace cis {
  template<typename NUM>
  bool between(NUM a, NUM b, NUM c) {
    if(a > c) swap(a,c);
    return a <= b && b <= c;
  }
}

bool in_seg(Seg l, P p) {
    bool in_x = cis::between(
	real(tail(l)), real(p),
        real(head(l)));
    bool in_y = cis::between(
	imag(tail(l)), imag(p),
        imag(head(l)));
    return in_x && in_y;
}

// general poly code -- not needed here
int ix_status;
// return true if segs intersect 
bool intersect(Seg a, Seg b, P& p) {
  P ab = pos(b) - pos(a);
  double r0 = cross(ab, dir(b));
  double s0 = cross(ab, dir(a));
  double d = cross(dir(a), dir(b));
  double r =  r0/d, s = s0/d;
  p = pos(a) + dir(a)*P(r);

  bool on_a = cis::between(0.,r,1.),
       on_b = cis::between(0.,s,1.);
  ix_status = 0; // 4 = par
  if(d == 0) ix_status |= 4;
  if(r0 == 0) ix_status |= 8;
  if(on_a) ix_status |= 1;
  if(on_b) ix_status |= 2;
  return on_a && on_b; // 8 = coll
}

int ccw(P a, P b, P c,
               bool coll = false) {
    P d = b-a, f = c-a;
    double diff = cross(d,f);
    if (diff < -e) return -1;
    if (diff >  e) return 1;
    if(coll) return 0;
    if( real(d)*real(f) < 0 ||
        imag(d)*imag(f) < 0) return -1;
    if(dot(d,d) < dot(f,f)) return 1;
    return 0;
}

bool onseg(P p, Seg l)
    { return 0 == ccw(tail(l), head(l), p); }

// consider segs_i([a,b], [c5,c2]) above:
// intersect() returns false (which is right).
//  Still, 2 of the 4 possible
// endpoint-on-other-seg will always be true,
//  so we can just test three of them
bool segs_intersect(Seg ab, Seg cd, P& i) {
    return intersect(ab,cd, i) ||
         ((ix_status & 12) == 12 &&
        (onseg(tail(cd), ab) ||
         onseg(head(cd), ab) ||
         onseg(tail(ab), cd)));
}
/// end infrastructure ////////

vector<P> towers;
vector<P> cities;

typedef int node;
typedef int weight;
typedef pair<node,node> edge;
typedef pair<edge,weight> nw;

vector<edge> roads; // vector of city pairs
    // map of city pairs to bool
map<edge, bool> roadmap;

    // map of city pairs to weight
map<edge, weight> roadweights;


P closest_point(Seg l, P p) {
    P x1 = tail(l);
    P r = p - x1;
    P q = x1 + proj(r, dir(l));
    // if q is off the line segment,
    // return the closest endpoint
    if(!in_seg(l, q)) {
        if(abs(q-head(l)) < abs(q-tail(l)))
            return head(l);
        return tail(l);
    }
    return q;
}

P closest_tower(P q) {
    double min_dist = 1000000;
    P closest = q;
    FOR(i,towers.size()) {
        double dist = abs(q - towers[i]);
        if(dist < min_dist) {
            min_dist = dist;
            closest = towers[i];
        }
    }
    return closest;
}

int B,C,R,Q;

int main() {
    double x,y;
  int m,n;
  int casenumber = 0;
  while(true) {
    // init data structures
    roads.clear();
    cities.clear();
    towers.clear();
    roadweights.clear();

    casenumber++;

    // read in lines
    cin >>B >> C >> R >> Q;

    if(B == 0 && C == 0
    && R == 0 && Q == 0) {
      return 0;
    }

    FOR(i,B) {
      cin >> x >> y;
      towers.push_back(P(x,y));
    }

    FOR(i,C) {
      cin >> x >> y;
      cities.push_back(P(x,y));
    }

    FOR(i,R) {
      cin >> m >> n;
      // decrease m and n by one,
      // to account for indexing
      --m; --n;
      roads.push_back(edge(m,n));
      roadmap[edge(m,n)] = true;
    }

    // compute weights for all roads
    // to do so, we'll take each line and
    // compute closest-point-on-lines (Q-points)
    // for each tower to the line.
    // Then, find closest tower for each Q-pt.
    // For a given line, the cardinality
    // of the set of closest towers,
    // minus one, is the weight of the line.
    FOR(road_no,R) {
      edge citypair = roads[road_no];
      Seg road = from_pts(cities[citypair.first],
                cities[citypair.second]);
      vector<P> qs;
      set<P> closest_towers = set<P>();
      FOR(tower_no, B) {
        P tower = towers[tower_no];
        P q = closest_point(road, tower);
        P ct = closest_tower(q);
        closest_towers.insert(ct);
      }
      int weight = closest_towers.size() - 1;
      // duplicate weights for city pairs
      roadweights[citypair] = weight;
      // no matter which city comes first
      roadweights[invert(citypair)] = weight;
    }

    map<edge, weight> spath;
    floyd(roads, roadweights, spath, C);

    cout << "Case " << casenumber << ":"<<endl;
    FOR(i,Q) {
      int s,d;
      cin >> s >> d;
      --s; --d;
      edge sd = edge(s,d);
      // max path can't be greater
      // than number of cities/nodes!
      if(spath[sd] > C) cout << "Impossible";
      else cout << spath[sd];
      cout << endl;
    }
  }
  return 0;
}
\end{lstlisting}
\normalsize
%%% Local Variables:
%%% mode: latex
%%% TeX-master: "ref-doc"
%%% End:
