C++ Piecewise::concat方法代码示例

    void draw(cairo_t *cr,
	      std::ostringstream *notify,
	      int width, int height, bool save) {
    
        D2<SBasis> B1 = b1_handle.asBezier();
        D2<SBasis> B2 = b2_handle.asBezier();
        Piecewise<D2<SBasis> >B;
        B.concat(Piecewise<D2<SBasis> >(B1));
        B.concat(Piecewise<D2<SBasis> >(B2));

// testing fuse_nearby_ends
        std::vector< Piecewise<D2<SBasis> > > pieces;
        pieces = fuse_nearby_ends(split_at_discontinuities(B),9);
        Piecewise<D2<SBasis> > C;
        for (unsigned i=0; i<pieces.size(); i++){
            C.concat(pieces[i]);
        }
// testing fuse_nearby_ends

        cairo_set_line_width (cr, .5);
        cairo_set_source_rgba (cr, 0., 0.5, 0., 1);
        //cairo_d2_sb(cr, B1);
        cairo_pw_d2_sb(cr, C);
        //cairo_pw_d2_sb(cr, B);
        cairo_stroke(cr);

        Piecewise<D2<SBasis> > uniform_B = arc_length_parametrization(B);
        cairo_set_source_rgba (cr, 0., 0., 0.9, 1);
        dot_plot(cr,uniform_B);
        cairo_stroke(cr);
        *notify << "pieces = " << uniform_B.size() << ";\n";

        Toy::draw(cr, notify, width, height, save);
    }        

/** Return a function which gives the angle of vect at each point.
 \param vect a piecewise parameteric curve.
 \param tol the maximum error allowed.
 \param order the maximum degree to use for approximation
 \relates Piecewise
*/
Piecewise<SBasis>
Geom::atan2(Piecewise<D2<SBasis> > const &vect, double tol, unsigned order){
    Piecewise<SBasis> result;
    Piecewise<D2<SBasis> > v = cutAtRoots(vect,tol);
    result.cuts.push_back(v.cuts.front());
    for (unsigned i=0; i<v.size(); i++){

        D2<SBasis> vi = RescaleForNonVanishingEnds(v.segs[i]);
        SBasis x=vi[0], y=vi[1];
        Piecewise<SBasis> angle;
        angle = divide (x*derivative(y)-y*derivative(x), x*x+y*y, tol, order);

        //TODO: I don't understand this - sign.
        angle = integral(-angle);
        Point vi0 = vi.at0(); 
        angle += -std::atan2(vi0[1],vi0[0]) - angle[0].at0();
        //TODO: deal with 2*pi jumps form one seg to the other...
        //TODO: not exact at t=1 because of the integral.
        //TODO: force continuity?

        angle.setDomain(Interval(v.cuts[i],v.cuts[i+1]));
        result.concat(angle);   
    }
    return result;
}

Geom::Piecewise<Geom::D2<Geom::SBasis> >
LPECopyRotate::doEffect_pwd2 (Geom::Piecewise<Geom::D2<Geom::SBasis> > const & pwd2_in)
{
    using namespace Geom;

    // I first suspected the minus sign to be a bug in 2geom but it is
    // likely due to SVG's choice of coordinate system orientation (max)
    start_pos = origin + dir * Rotate(-deg_to_rad(starting_angle)) * dist_angle_handle;
    double rotation_angle_end = rotation_angle;
    if(copiesTo360){
        rotation_angle_end = 360.0/(double)num_copies;
    }
    rot_pos = origin + dir * Rotate(-deg_to_rad(starting_angle + rotation_angle_end)) * dist_angle_handle;

    A = pwd2_in.firstValue();
    B = pwd2_in.lastValue();
    dir = unit_vector(B - A);

    Piecewise<D2<SBasis> > output;

    Affine pre = Translate(-origin) * Rotate(-deg_to_rad(starting_angle));
    for (int i = 0; i < num_copies; ++i) {
        // I first suspected the minus sign to be a bug in 2geom but it is
        // likely due to SVG's choice of coordinate system orientation (max)
        Rotate rot(-deg_to_rad(rotation_angle_end * i));
        Affine t = pre * rot * Translate(origin);
        output.concat(pwd2_in * t);
    }

    return output;
}

/** Compute the cosine of a function.
 \param f function
 \param tol maximum error
 \param order maximum degree polynomial to use
*/
Piecewise<SBasis> cos(Piecewise<SBasis> const &f, double tol, int order){
    Piecewise<SBasis> result;
    for (unsigned i=0; i<f.size(); i++){
        Piecewise<SBasis> cosfi = cos(f.segs[i],tol,order);
        cosfi.setDomain(Interval(f.cuts[i],f.cuts[i+1]));
        result.concat(cosfi);
    }
    return result;
}

/** Reparameterise M to have unit speed.
 \param M the Element.
 \param tol the maximum error allowed.
 \param order the maximum degree to use for approximation
 \relates Piecewise
*/
Piecewise<D2<SBasis> >
Geom::arc_length_parametrization(Piecewise<D2<SBasis> > const &M,
                                 unsigned order,
                                 double tol){
    Piecewise<D2<SBasis> > result;
    for (unsigned i=0; i<M.size(); i++) {
        result.concat( arc_length_parametrization(M[i],order,tol) );
    }
    return result;
}

/** Compute the sqrt of a function.
 \param f function
*/
Piecewise<SBasis> sqrt(Piecewise<SBasis> const &f, double tol, int order){
    Piecewise<SBasis> result;
    Piecewise<SBasis> zero = Piecewise<SBasis>(Linear(tol*tol));
    zero.setDomain(f.domain());
    Piecewise<SBasis> ff=max(f,zero);

    for (unsigned i=0; i<ff.size(); i++){
        Piecewise<SBasis> sqrtfi = sqrt_internal(ff.segs[i],tol,order);
        sqrtfi.setDomain(Interval(ff.cuts[i],ff.cuts[i+1]));
        result.concat(sqrtfi);
    }
    return result;
}

/** Return a Piecewise<D2<SBasis> > which points in the same direction as V_in, but has unit_length.
 \param V_in the original path.
 \param tol the maximum error allowed.
 \param order the maximum degree to use for approximation

unitVector(x,y) is computed as (b,-a) where a and b are solutions of:
     ax+by=0 (eqn1)   and   a^2+b^2=1 (eqn2)

 \relates Piecewise
*/
Piecewise<D2<SBasis> >
Geom::unitVector(Piecewise<D2<SBasis> > const &V, double tol, unsigned order){
    Piecewise<D2<SBasis> > result;
    Piecewise<D2<SBasis> > VV = cutAtRoots(V);
    result.cuts.push_back(VV.cuts.front());
    for (unsigned i=0; i<VV.size(); i++){
        Piecewise<D2<SBasis> > unit_seg;
        unit_seg = unitVector(VV.segs[i],tol, order);
        unit_seg.setDomain(Interval(VV.cuts[i],VV.cuts[i+1]));
        result.concat(unit_seg);   
    }
    return result;
}

Geom::Piecewise<Geom::D2<Geom::SBasis> >
doEffect_pwd2 (Geom::Piecewise<Geom::D2<Geom::SBasis> > & pwd2_in, Geom::Piecewise<Geom::D2<Geom::SBasis> > & pattern)
{
    using namespace Geom;

    Piecewise<D2<SBasis> > uskeleton = arc_length_parametrization(pwd2_in, 2, .1);
    uskeleton = remove_short_cuts(uskeleton,.01);
    Piecewise<D2<SBasis> > n = rot90(derivative(uskeleton));
    n = force_continuity(remove_short_cuts(n,.1));

    D2<Piecewise<SBasis> > patternd2 = make_cuts_independent(pattern);
    Piecewise<SBasis> x = Piecewise<SBasis>(patternd2[0]);
    Piecewise<SBasis> y = Piecewise<SBasis>(patternd2[1]);
    Interval pattBnds = *bounds_exact(x);
    x -= pattBnds.min();
    Interval pattBndsY = *bounds_exact(y);
    y -= (pattBndsY.max()+pattBndsY.min())/2;

    int nbCopies = int(uskeleton.cuts.back()/pattBnds.extent());
    double scaling = 1;

    double pattWidth = pattBnds.extent() * scaling;

    if (scaling != 1.0) {
        x*=scaling;
    }

    double offs = 0;
    Piecewise<D2<SBasis> > output;
    for (int i=0; i<nbCopies; i++){
        output.concat(compose(uskeleton,x+offs)+y*compose(n,x+offs));
        offs+=pattWidth;
    }

    return output;
}

// Main effect body...
Geom::Piecewise<Geom::D2<Geom::SBasis> >
LPESketch::doEffect_pwd2 (Geom::Piecewise<Geom::D2<Geom::SBasis> > const & pwd2_in)
{
    using namespace Geom;
    //If the input path is empty, do nothing.
    //Note: this happens when duplicating a 3d box... dunno why.
    if (pwd2_in.size()==0) return pwd2_in;

    Piecewise<D2<SBasis> > output;

    // some variables for futur use (for construction lines; compute arclength only once...)
    // notations will be : t = path time, s = distance from start along the path.
    Piecewise<SBasis> pathlength;
    double total_length = 0;

    //TODO: split Construction Lines/Approximated Strokes into two separate effects?

    //----- Approximated Strokes.
    std::vector<Piecewise<D2<SBasis> > > pieces_in = split_at_discontinuities (pwd2_in);
    
    //work separately on each component.
    for (unsigned pieceidx = 0; pieceidx < pieces_in.size(); pieceidx++){

        Piecewise<D2<SBasis> > piece = pieces_in[pieceidx];
        Piecewise<SBasis> piecelength = arcLengthSb(piece,.1);
        double piece_total_length = piecelength.segs.back().at1()-piecelength.segs.front().at0();
        pathlength.concat(piecelength + total_length);
        total_length += piece_total_length;
        

        //TODO: better check this on the Geom::Path.
        bool closed = piece.segs.front().at0() == piece.segs.back().at1(); 
        if (closed){ 
            piece.concat(piece);
            piecelength.concat(piecelength+piece_total_length);
        }

        for (unsigned i = 0; i<nbiter_approxstrokes; i++){
            //Basic steps: 
            //- Choose a rdm seg [s0,s1], find coresponding [t0,t1], 
            //- Pick a rdm perturbation delta(s), collect 'piece(t)+delta(s(t))' over [t0,t1] into output.

            // pick a point where to start the stroke (s0 = dist from start).
            double s1=0.,s0 = ends_tolerance*strokelength+0.0001;//the root finder might miss 0.  
            double t1, t0;
            double s0_initial = s0;
            bool done = false;// was the end of the component reached?

            while (!done){
                // if the start point is already too far... do nothing. (this should not happen!)
                if (!closed && s1>piece_total_length - ends_tolerance.get_value()*strokelength) break;
                if ( closed && s0>piece_total_length + s0_initial) break;

                std::vector<double> times;  
                times = roots(piecelength-s0);  
                t0 = times.at(0);//there should be one and only one solution!!
                
                // pick a new end point (s1 = s0 + strokelength).
                s1 = s0 + strokelength*(1-strokelength_rdm);
                // don't let it go beyond the end of the orgiginal path.
                // TODO/FIXME: this might result in short strokes near the end...
                if (!closed && s1>piece_total_length-ends_tolerance.get_value()*strokelength){
                    done = true;
                    //!!the root solver might miss s1==piece_total_length...
                    if (s1>piece_total_length){s1 = piece_total_length - ends_tolerance*strokelength-0.0001;}
                }
                if (closed && s1>piece_total_length + s0_initial){
                    done = true;
                    if (closed && s1>2*piece_total_length){
                        s1 = 2*piece_total_length - strokeoverlap*(1-strokeoverlap_rdm)*strokelength-0.0001;
                    }
                }
                times = roots(piecelength-s1);  
                if (times.size()==0) break;//we should not be there.
                t1 = times[0];
                
                //pick a rdm perturbation, and collect the perturbed piece into output.
                Piecewise<D2<SBasis> > pwperturb = computePerturbation(s0-0.01,s1+0.01);
                pwperturb = compose(pwperturb,portion(piecelength,t0,t1));

                output.concat(portion(piece,t0,t1)+pwperturb);

                //step points: s0 = s1 - overlap.
                //TODO: make sure this has to end?
                s0 = s1 - strokeoverlap*(1-strokeoverlap_rdm)*(s1-s0);
            }
        }
    }

    //----- Construction lines.
    //TODO: choose places according to curvature?.

    //at this point we should have:
    //pathlength = arcLengthSb(pwd2_in,.1);
    //total_length = pathlength.segs.back().at1()-pathlength.segs.front().at0();
    Piecewise<D2<SBasis> > m = pwd2_in;
    Piecewise<D2<SBasis> > v = derivative(pwd2_in);
    Piecewise<D2<SBasis> > a = derivative(v);

//.........这里部分代码省略.........

    void draw(cairo_t *cr,
              std::ostringstream *notify,
              int width, int height, bool save, std::ostringstream *timer_stream) {

        if (first_time)
        {
            first_time = false;
            sliders[0].geometry(Point(50, 50), 100);
        }
        size_t const num_points = static_cast<size_t>(sliders[0].value());

        D2<SBasis> B1 = b_handle.asBezier();
        Piecewise<D2<SBasis> >B;
        B.concat(Piecewise<D2<SBasis> >(B1));

        // testing fuse_nearby_ends
        std::vector< Piecewise<D2<SBasis> > > pieces;
        pieces = fuse_nearby_ends(split_at_discontinuities(B),9);
        Piecewise<D2<SBasis> > C;
        for (unsigned i=0; i<pieces.size(); i++){
            C.concat(pieces[i]);
        }
        // testing fuse_nearby_ends

        cairo_set_line_width (cr, 2.);
        cairo_set_source_rgba (cr, 0., 0.5, 0., 1);
        //cairo_d2_sb(cr, B1);
        //cairo_pw_d2_sb(cr, C);
        //cairo_pw_d2_sb(cr, B);
        cairo_stroke(cr);

        Timer tm;
        Timer::Time als_time = tm.lap();

        cairo_set_source_rgba (cr, 0., 0., 0.9, 1);
        //dot_plot(cr,uniform_B);
        cairo_stroke(cr);

        std::cout << B[0] << std::endl;

        Geom::Affine translation;

        Geom::Path original_path;
        //original_bezier.append(B[0]);
        //original_bezier.appendNew<CubicBezier> (B[0]);
        CubicBezier original_bezier(b_handle.pts);
        original_path.append(original_bezier);

        std::vector<double> initial_t;
        std::vector<Geom::Point> curve_points;
        if (randomize_times) {
            std::uniform_real_distribution<double> dist_t(0,1);
            for (size_t ii = 0; ii < num_points; ++ii) {
                double const t = dist_t(generator);
                initial_t.push_back(t);
            }
            std::sort(initial_t.begin(), initial_t.end());
            double const min = initial_t.front();
            double const max = initial_t.back();
            for (auto& t : initial_t) {
                t = (t-min)/(max-min);
            }
            for (auto const t : initial_t) {
                curve_points.push_back(original_bezier.pointAt(t));
            }
        }
        else {
            for (size_t ii = 0; ii < num_points; ++ii) {
                double const t = static_cast<double>(ii) / (num_points-1);
                Geom::Point const p = original_bezier.pointAt(t);
                initial_t.push_back(t);
                curve_points.push_back(p);
            }
        }

        cairo_set_source_rgba (cr, 0., 0., .9, 1);
        cairo_path(cr, original_path);
        draw_text(cr, original_path.initialPoint(), "original curve and old fit");


        Geom::CubicBezier fitted_new;
        Geom::CubicBezier fitted_new_a;
        Geom::Point very_old_version_raw[4];
        bezier_fit_cubic(very_old_version_raw, curve_points.data(), curve_points.size(), 0.);
        Geom::CubicBezier very_old_bezier(
                    very_old_version_raw[0],
                very_old_version_raw[1],
                very_old_version_raw[2],
                very_old_version_raw[3]
                );

        Geom::Path very_old_version_path;
        very_old_version_path.append(very_old_bezier);

        cairo_set_source_rgba (cr, .7, .7, 0., 1);
        cairo_stroke(cr);
        cairo_path(cr, very_old_version_path);

        cairo_set_source_rgba (cr, 0., 0., .9, 1);
        cairo_stroke(cr);
//.........这里部分代码省略.........

Geom::Piecewise<Geom::D2<Geom::SBasis> >
LPERecursiveSkeleton::doEffect_pwd2 (Geom::Piecewise<Geom::D2<Geom::SBasis> > const & pwd2_in)
{
    using namespace Geom;

    Piecewise<D2<SBasis> > output;
    std::vector<Piecewise<D2<SBasis> > > pre_output;

    double prop_scale = 1.0;

    D2<Piecewise<SBasis> > patternd2 = make_cuts_independent(pwd2_in);
    Piecewise<SBasis> x0 = false /*vertical_pattern.get_value()*/ ? Piecewise<SBasis>(patternd2[1]) : Piecewise<SBasis>(patternd2[0]);
    Piecewise<SBasis> y0 = false /*vertical_pattern.get_value()*/ ? Piecewise<SBasis>(patternd2[0]) : Piecewise<SBasis>(patternd2[1]);
    OptInterval pattBndsX = bounds_exact(x0);
    OptInterval pattBndsY = bounds_exact(y0);

    if ( !pattBndsX || !pattBndsY) {
        return pwd2_in;
    }

    x0 -= pattBndsX->min();
    y0 -= pattBndsY->middle();

    double xspace  = 0;//spacing;
    double noffset = 0;//normal_offset;
    double toffset = 0;//tang_offset;
    if (false /*prop_units.get_value()*/){
        xspace  *= pattBndsX->extent();
        noffset *= pattBndsY->extent();
        toffset *= pattBndsX->extent();
    }

    y0+=noffset;

    output = pwd2_in;

    for (int i = 0; i < iterations; ++i) {
        std::vector<Piecewise<D2<SBasis> > > skeleton = split_at_discontinuities(output);

        output.clear();
        for (unsigned idx = 0; idx < skeleton.size(); idx++){
            Piecewise<D2<SBasis> > path_i = skeleton[idx];
            Piecewise<SBasis> x = x0;
            Piecewise<SBasis> y = y0;
            Piecewise<D2<SBasis> > uskeleton = arc_length_parametrization(path_i,2,.1);
            uskeleton = remove_short_cuts(uskeleton,.01);
            Piecewise<D2<SBasis> > n = rot90(derivative(uskeleton));
            n = force_continuity(remove_short_cuts(n,.1));

            double scaling = 1;
            scaling = (uskeleton.domain().extent() - toffset)/pattBndsX->extent();

            // TODO investigate why pattWidth is not being used:
            double pattWidth = pattBndsX->extent() * scaling;

            if (scaling != 1.0) {
                x*=scaling;
            }

            if ( true /*scale_y_rel.get_value()*/ ) {
                y*=(scaling*prop_scale);
            } else {
                if (prop_scale != 1.0) y *= prop_scale;
            }
            x += toffset;

            output.concat(compose(uskeleton,x)+y*compose(n,x));
        }
    }

    return output;
}