Problem: E --Disappear Drive

Story

The person in D likes basketball, but he likes shooting, so other techniques are crazy. I'm not particularly good at dribbling, and when I enter the opponent's defensive range, I always get the ball stolen. So I decided to come up with a deadly dribble that would surely pull out any opponent. After some effort, he finally completed the disappearing drive, "Disapia Drive".

This Apia Drive creates a gap by provoking the opponent to deprive him of his concentration and take his eyes off, and in that gap he goes beyond the wall of dimension and slips through the opponent's defensive range. It's a drive. However, crossing the dimensional wall puts a heavy burden on the body, so there is a limit to the number of times the dimensional wall can be crossed. Also, because he always uses his concentration, he can only change direction once, including normal movements. How do you move on the court to reach the goal in the shortest time without being robbed of the ball?

Problem

Consider a rectangle on a two-dimensional plane with (0,0) at the bottom left and (50,94) at the top right. In addition, there are N circles on this plane, the center position of the i-th circle is (x_i, y_i), and the radius is r_i. There is no overlap between the circles.

Let (25,0) be the point S and (25,94) be the point G, and consider the route from S to G. A "path" consists of (1) a line segment connecting S and G, or (2) a combination of a line segment SP and a line segment GP at any point P inside a rectangle. In the middle of the route from S to G, the period from entering the inside of a circle to exiting the inside of the circle is defined as "the section passing through the inside of the circle". However, it is assumed that the circumferences of the rectangle and the circle are not included in the rectangle and the circle, respectively.

Find the length of the shortest route among the routes where the number of sections passing through the inside of the circle is D or less. If you can't reach the goal, return -1.

Input

The input consists of the following format.


N D
x_1 y_1 r_1
...
x_N y_N r_N

The first line consists of two integers, and the number N of circles and the number of times D that can pass inside the circle are given, separated by one blank character. The following N lines are given circle information. The i + 1 line (1 \ leq i \ leq N) consists of three integers, and the x-coordinate x_i, y-coordinate y_i, and radius r_i of the center of the i-th circle are given by separating them with one blank character.

Constraints:

* 0 \ leq N \ leq 5
* 0 \ leq D \ leq 5
* 0 \ leq x_i \ leq 50
* 0 \ leq y_i \ leq 94
* 1 \ leq r_i \ leq 100
* It can be assumed that any two circles are separated by 10 ^ {-5} or more.
* It can be assumed that S and G are not included inside any circle.
* It can be assumed that S and G are separated from any circle by 10 ^ {-5} or more.
* It can be assumed that the point P in the shortest path is more than 10 ^ {-5} away from the circumference of the rectangle and any circumference.



Output

Output the length of the shortest path on one line. Be sure to start a new line at the end of the line. However, if you cannot reach the goal, return -1. Absolute error or relative error of 10 ^ {-7} or less is allowed for the correct answer.

Sample Input 1


Ten
25 47 10

Sample Output 1


96.2027355887

DisappearDrive_sample1.png

The detour route is the shortest because it cannot disappear.

Sample Input 2


1 1
25 47 10

Sample Output 2


94.0000000000

DisappearDrive_sample2.png

It can disappear, so you just have to go straight through it.

Sample Input 3


Ten
20 47 5

Sample Output 3


94.0000000000

DisappearDrive_sample3.png

Since the circumference is not included in the circle, it can pass through without disappearing.

Sample Input 4


Ten
25 47 40

Sample Output 4


-1

DisappearDrive_sample4.png

You can't reach the goal because you can't disappear or detour.

Sample Input 5


5 2
11 10 16
33 40 18
20 66 10
45 79 14
22 85 8

Sample Output 5


96.1320937224

DisappearDrive_sample5.png





Example

Input

1 0
25 47 10


Output

96.2027355887

Solution (in CPP):

Solution:
#include <bits/stdc++.h>
#define REP(i,n) for(int i=0; i<(int)(n); ++i)

using namespace std;

using namespace std;

typedef complex<double> P;
const double EPS = 1e-8;

// 誤差を加味した符号判定
int sign(double a, double b = 0){
    if(a - b > EPS) return +1;
    if(a - b < -EPS) return -1;
    return 0;
}

// 比較演算子
namespace std{
    bool operator < (const P& a, const P& b) {
        return real(a) != real(b) ? real(a) < real(b) : imag(a) < imag(b);
    }
}

// 内積・外積
double dot(P a, P b){
    return real(conj(a) * b);
}
double cross(P a, P b){
    return imag(conj(a) * b);
}

// OAとOBのなす符号付き角度 [-pi, pi]
// example : (1, 0), (0, 1) -> pi/2
double angle(P a, P b){
    return arg(conj(a) * b);
}

// aをc中心にb[rad]回転
// verify : not yet.
P rotate(P a, double b, P c = P()){
    return (a - c) * polar(1.0, b) + c;
}

// 直線ABに対する点Cの位置
int ccw(P a, P b, P c) {
    b -= a; c -= a;
    if (cross(b, c) > +EPS)   return +1; // 反時計回り
    if (cross(b, c) < -EPS)   return -1; // 時計回り
    if (dot(b, c) < 0)     return +2; // c--a--b の順番で一直線上
    if (norm(b) < norm(c)) return -2; // a--b--c の順番で一直線上
    return 0;                         // 点が線分ab上にある
}

enum{ OUT, ON, IN };
struct L : public vector<P> {
    L(){}
    L(P a, P b) {
        push_back(a); push_back(b);
    }
    P vector() const {
        return back() - front();
    }
};

// 注意: 端点で交わったり直線が重なったりする場合も交差していると判定する

// 二直線の平行判定
// verify : aoj0021
bool paralell(L l, L m){
    return sign(cross(l.vector(), m.vector())) == 0;
}

// 二直線の同一判定
bool equalLL(L l, L m){
    return sign(cross(l.vector(), m[0] - l[0])) == 0;
}

// 直線と点の交差判定
bool iLP(L l, P p) {
    // 直線lとl[0]からpへの直線が平行
    return sign(cross(l.vector(), p - l[0])) == 0;
}

// 線分と点の交差判定(端点の処理に注意)
// verify : aoj1279
bool iSP(L s, P p) {
    return ccw(s[0], s[1], p) == 0;
}

// 直線と線分の交差判定(線分が重なっている時に注意)
bool iLS(L l, L s) {
    // 直線lについて、線分sの端点が異なる側にある
    return sign(cross(l.vector(), s[0] - l[0]) * cross(l.vector(), s[1] - l[0])) <= 0;
}

// 二つの線分の交差判定(線分が重なっている時や端点の処理に注意)
bool iSS(L s, L t) {
    return ccw(s[0], s[1], t[0]) * ccw(s[0], s[1], t[1]) <= 0 &&
        ccw(t[0], t[1], s[0]) * ccw(t[0], t[1], s[1]) <= 0;
}


// 点pから直線lに対する射影
P proj(L l, P p){
    double t = dot(p - l[0], l.vector()) / norm(l.vector());
    return l[0] + t * l.vector();
}

// 点pの直線lに関する反射
P refl(L l, P p){
    return 2.0 * proj(l, p) - p;
}

// 直線と点の距離
double dLP(L l, P p){
    // return abs(p - projection(l, p));
    return abs(cross(l.vector(), p - l[0])) / abs(l.vector());
}

// 線分と点の距離 ( not verified !!! )
double dSP(L s, P p){
    if(sign(dot(s.vector(), p - s[0])) <= 0) return abs(p - s[0]);
    if(sign(dot(-s.vector(), p - s[1])) <= 0) return abs(p - s[1]);
    return dLP(s, p);
}

// 直線と直線の距離
double dLL(L l, L m){
    // 平行でないときは0, 平行のときは垂線の長さ
    return paralell(l, m) ? dLP(l, m[0]) : 0;
}

// 直線と線分の距離
double dLS(L l, L s){
    if(iLS(l, s)) return 0;
    return min(dLP(l, s[0]), dLP(l, s[1]));
}

// 線分と線分の距離
double dSS(L s, L t){
    if(iSS(s, t)) return 0;
    return min({dSP(s, t[0]), dSP(s, t[1]), dSP(t, s[0]), dSP(t, s[1])});
}

// 直線と直線の交点
P pLL(L l, L m){
    double A = cross(l.vector(), m.vector());
    double B = cross(l.vector(), l[1] - m[0]);
    if(sign(A) == 0 && sign(B) == 0) return m[0]; // 二直線が重なっている
    if(sign(A) == 0) assert(false); // 直線が交わらない
    return m[0] + m.vector() * B / A;
}

typedef vector<P> Pol; // 反時計回りを仮定

// Polの要素へのアクセス
P curr(const Pol& a, int x){ return a[x]; }
P next(const Pol& a, int x){ return a[(x + 1) % a.size()]; }
P prev(const Pol& a, int x){ return a[(x - 1 + a.size()) % a.size()]; }

// 点が多角形のどこにあるのか判定する
// verify : aoj0012
int contains(const Pol& A, P p){
    // 点pから半直線をひき、辺と交差する回数を数える
    bool in = false;
    for(int i = 0; i < A.size(); i++){
        P a = curr(A, i) - p;
        P b = next(A, i) - p;
        if(a.imag() > b.imag()) swap(a, b);
        // aからbの直線がy=0と交わり、その交点は原点の右側である
        if(a.imag() <= 0 && 0 < b.imag() && cross(a, b) < 0){ 
            in = !in;
        }
        if(sign(cross(a, b)) == 0 && sign(dot(a, b)) <= 0) return ON;
    }
    return in ? IN : OUT;
}

// 多角形の面積
// verify : aoj0079 aoj1100
double area(const Pol& A) {
    double res = 0;
    for(int i = 0; i < A.size(); i++){
        res += cross(curr(A, i), next(A, i));
    }
    return abs(res) / 2.0;
}

// 凸包
Pol convex_hull(vector<P> ps) {
    int n = ps.size(), k = 0;
    sort(ps.begin(), ps.end());
    vector<P> ch(2*n);
    for (int i = 0; i < n; ch[k++] = ps[i++]){ // lower-hull
        while (k >= 2 && ccw(ch[k-2], ch[k-1], ps[i]) <= 0) --k;
    }
    for (int i = n-2, t = k+1; i >= 0; ch[k++] = ps[i--]){ // upper-hull
        while (k >= t && ccw(ch[k-2], ch[k-1], ps[i]) <= 0) --k;
    }
    ch.resize(k-1);
    return ch;
}

bool is_convex(const Pol& A){
    for(int i = 0; i < A.size(); i++){
        if(ccw(prev(A, i), curr(A, i), next(A, i)) > 0) return false;
    }
    return true;
}

// 凸多角形の直線による切断。直線の左側だけ残す
// verify : aoj1283
Pol convex_cut(const Pol& A, L l){
    Pol B;
    for(int i = 0; i < A.size(); i++){
        P a = curr(A, i), b = next(A, i);
        if(ccw(l[0], l[1], a) != -1) B.push_back(a); //Aが直線lの右側でない
        if(ccw(l[0], l[1], a) * ccw(l[0], l[1], b) < 0)
            B.push_back(pLL(l, L(a, b)));
    }
    return B;
}
// 垂直二等分線
// verify: maximamcup2013 D
L bisector(P a, P b){
    P mid = (a + b) / 2.0;
    P vec = (mid - a) * P(0.0, 1.0);
    return L(mid, mid + vec);
}
// 点集合psのうちs番目のボロノイ領域
// verify: maximamcup2013 D
Pol voronoi_cell(Pol A, const vector<P>& ps, int s){
    for(int i = 0; i < ps.size(); i++){
        if(i != s) A = convex_cut(A, bisector(ps[s], ps[i]));
    }
    return A;
}
struct C {
    P p;
    double r;
    C() {}
    C(P p, double r) : p(p), r(r) { }
};

// 円と点の内外判定
int contains(C c, P p){
    double d = abs(c.p - p);
    if(sign(d - c.r) > 0) return OUT;
    if(sign(d - c.r) == 0) return ON;
    return IN;
}

// 円と線分の交差判定(境界を含む)
// Verified: AOJ 0129
bool iCS(C c, L l){
    int c1 = contains(c, l[0]);
    int c2 = contains(c, l[1]);
    if(c1 > c2) swap(c1, c2);

    // (OUT, OUT) (OUT, ON) (OUT, IN) (ON, ON) (ON, IN) (IN, IN) の6通り
    if(c1 == OUT && c2 == IN) return true;
    if(c1 == IN  && c2 == IN) return false;
    if(c1 == ON) return true; // (接するとき) 
    double d = dSP(l, c.p);
    if(sign(d - c.r) < 0) return true;
    if(sign(d - c.r) == 0) return true; // (接するとき)
    if(sign(d - c.r) > 0) return false;
    assert(false);
}

// 二つの円の交差判定(接する時を含む)
bool iCC(C c, C d){
    // 円の中心同士の距離が、半径の和以下であり、半径の差以上である
    double e = abs(d.p - d.p);
    return sign(e - (c.r + d.r)) <= 0 && sign(e - abs(c.r - d.r)) >= 0;
}

// 円と直線の交点
// verify : aoj2045
vector<P> pLC(L l, C c) {
    vector<P> res;
    P center = proj(l, c.p);
    double d = abs(center - c.p);
    double tt = c.r * c.r - d * d;
    if(tt < 0 && tt > -EPS) tt = 0;
    if(tt < 0) return res;
    double t = sqrt(tt);
    P vect = l.vector();
    vect /= abs(vect);
    res.push_back(center - vect * t);
    if (t > EPS) {
        res.push_back(center + vect * t);
    }
    return res;
}

// 円と線分の交点
vector<P> pSC(L s, C c) {
    vector<P> ret;
    vector<P> nret = pLC(s, c);
    for (int i = 0; i < nret.size(); i++) {
        if (iSP(s, nret[i])) ret.push_back(nret[i]);
    }
    return ret;
}

// 円と円の交点
// verify : aoj1183
vector<P> pCC(C a, C b){
    vector<P> res;

    double l = abs(b.p - a.p);

    if(sign(l) == 0 && sign(a.r - b.r) == 0) assert(false); // 解が無限に存在する
    if(sign(l - abs(a.r - b.r)) < 0 || sign(l - (a.r + b.r)) > 0) return res; // 解が存在しない

    double th1 = arg(b.p - a.p);
    if(sign(l - abs(a.r - b.r)) == 0 || sign(l - (a.r + b.r)) == 0){
        res.push_back(a.p + polar(a.r, th1));
    }else {
        double th2 = acos( (a.r * a.r - b.r * b.r + l * l) / (2 * a.r * l) );
        res.push_back(a.p + polar(a.r, th1 - th2));
        res.push_back(a.p + polar(a.r, th1 + th2));
    }
    return res;
}

// 2点を通る半径rの円の中心
// verify : aoj1132
vector<P> touching_circle2(P a, P b, double r){
    vector<P> res;

    double d = abs(b - a);
    if(d > 2 * r) return res;

    P mid = 0.5 * (a + b);
    P dir = polar(sqrt(r * r - d * d / 4), arg(b - a) + M_PI / 2);
    res.push_back(mid + dir);
    res.push_back(mid - dir);
    return res;
}

// 3点を通る円
C touching_circle3(P a, P b, P c){
    // ２つの垂直二等分線の交点が円の中心
    P mid_ab = (a + b) / 2.0;
    L bis_ab(mid_ab, (mid_ab - a) * P(0.0, 1.0));
    P mid_bc = (b + c) / 2.0;
    L bis_bc(mid_bc, (mid_bc - b) * P(0.0, 1.0));

    assert(!paralell(bis_ab, bis_bc)); 

    P center = pLL(bis_ab, bis_bc);
    return C(center, abs(a - center));
}

// 円と円の共通部分の面積を求める．
// ref: nya3j
double cc_area(C c1, C c2) {
    double d = abs(c1.p - c2.p);
    if (c1.r + c2.r < d + EPS) {
        return 0.0;
    } else if (d < abs(c1.r - c2.r) + EPS) {
        double r = min(c1.r, c2.r); // 元は c1.r >? c2.r だった．
        return r * r * M_PI;
    } else {
        double rc = (d*d + c1.r*c1.r - c2.r*c2.r) / (2*d);
        double theta = acos(rc / c1.r);
        double phi = acos((d - rc) / c2.r);
        return c1.r*c1.r*theta + c2.r*c2.r*phi - d*c1.r*sin(theta);
    }
}
// 円の接線 (中心から偏角thの点で接する接線)
// verified: AOJ 2201 Immortal Jewels 
L circle_tangent(const C& c, double th){
    P p0 = c.p + polar(c.r, th);
    P p1 = p0 + polar(1.0, th + M_PI / 2);
    return L(p0, p1);
}

// 二つの円の共通接線 (cの中心から接点へのベクトルの偏角を返す)
// verified: AOJ 2201 Immortal Jewels 
// 参考: http://geom.web.fc2.com/geometry/circle-circle-tangent.html
vector<double> common_tangents(const C& c, const C& d){
    vector<double> res;
    P v = d.p - c.p;
    double l = abs(v); // 二円の中心間の距離
    double a = arg(v); // 二円の中心間の偏角
    if(sign(l - abs(c.r - d.r)) > 0){
        // 交わる or 外接 or 離れている
        // 二つの外側接線
        double a1 = acos((c.r - d.r) / l);
        res.push_back(a + a1);
        res.push_back(a - a1);
        if(sign(l - (c.r + d.r)) > 0){
            // 離れている
            // 二つの内側接線
            double a2 = acos((c.r + d.r) / l);
            res.push_back(a + a2);
            res.push_back(a - a2);
        }
    }
    if((sign(l - abs(c.r - d.r)) == 0 || sign(l - (c.r + d.r)) == 0) && sign(l) != 0){
        // 内接 or 外接
        // 一つの接線
        res.push_back(a);
    }
    return res;
}

// 1点を通る円の接線( pがcの外側にあることが前提条件 )
// verified : AOJ 2579
vector<L> tangents_through_point(const C& c, const P& p){
    vector<L> tangents;
    double d = abs(c.p - p);
    // d ^ 2 == r ^ 2 + e ^ 2
    double e = sqrt(d * d - c.r * c.r); // 点pと円の接点の距離
    // d * sin(th) = r
    double th = asin(c.r / d);
    P q1 = p + (c.p - p) * polar(1.0, +th) * e / d;
    P q2 = p + (c.p - p) * polar(1.0, -th) * e / d;
    tangents.push_back(L(p, q1));
    tangents.push_back(L(p, q2));
    return tangents;
}
// 三角形の内心
// verify : aoj 1301
P incenter(P p1, P p2, P p3){
    double a = abs(p2 - p3);
    double b = abs(p3 - p1);
    double c = abs(p1 - p2);
    return (a * p1 + b * p2 + c * p3) / (a + b + c);
}

int main(){
    const double W = 50;
    const double H = 94;
    int N, D;
    cin >> N >> D;
    vector<C> c(N);
    REP(i, N){
        double x, y, r;
        cin >> x >> y >> r;
        c[i] = {{x, y}, r};
    }
    const P start(W / 2, 0);
    const P goal (W / 2, H);
    vector<L> st, gt;
    for(int i = 0; i < N; i++){
        vector<L> sv = tangents_through_point(c[i], start);
        vector<L> gv = tangents_through_point(c[i], goal );
        st.insert(st.end(), sv.begin(), sv.end());
        gt.insert(gt.end(), gv.begin(), gv.end());
    }
    const double INF = 1e9;
    double ans = INF;

    // 直線で進んだとき
    {
        L shortest(start, goal);
        int cnt = 0;
        for(int i = 0; i < N; i++){
            if(iCS(c[i], shortest)) cnt++;
        }
        if(cnt <= D){
            ans = abs(shortest.vector());
        }
    }

    // 円を接するときの経路
    for(L sl : st){
        for(L gl : gt){
            if(paralell(sl, gl)) continue;
            P cp = pLL(sl, gl);
            if(!(sign(0, cp.real()) <= 0 && sign(cp.real() ,W) <= 0)) continue;
            if(!(sign(0, cp.imag()) <= 0 && sign(cp.imag(), H) <= 0)) continue;
            L sm(start, cp);
            L gm(cp, goal);
            // cout << sm[0] << " -> " << sm[1] << endl;
            // cout << gm[0] << " -> " << gm[1] << endl;
            int cnt = 0;
            for(int i = 0; i < N; i++){
                if(iCS(c[i], sm) && sign(dSP(sm, c[i].p), c[i].r) != 0) cnt++;
                if(iCS(c[i], gm) && sign(dSP(gm, c[i].p), c[i].r) != 0) cnt++;
                if(contains(c[i], cp) == IN) cnt--;
            }
            if(cnt <= D){
                double d = abs(sm.vector()) + abs(gm.vector());
                ans = min(ans, d);
            }
        }
    }

    if(ans == INF) {
        printf("-1\n");
    } else {
        printf("%.12f\n", ans);
    }
    return 0;
}