ivp_mindist_debug.cxx

hl2 source code. Do not use it illegal.

    m_cache_0->remove_reference();
    m_cache_1->remove_reference();
    
    // compare brute force result with original mindist state
    IVP_DOUBLE dist_0 = sqrt(min_qdist) - mindist->sum_extra_radius;
    IVP_DOUBLE dist_1 = mindist->get_length();
    IVP_DOUBLE diff = IVP_Inline_Math::fabsd(dist_0 - dist_1) + P_DOUBLE_EPS;
    
    IVP_DOUBLE perc = diff * 100.0f / dist_1;
    // check if difference is neglectible
    if(perc < 50.0f){
	return IVP_OK;
    }
    
    // somewhat greater differences ...
    printf("\nBrute Force Recalculation of Mindist shows differences:\n");
    printf("Dsts: md %f, bf %f, diff: %f percent: %f\n", dist_0, dist_1, diff, perc);
    printf("Coords of synapses (original 0/1, brute force 0/1):\n");
    
    //IVP_DOUBLE qd0 = IVP_CLS.calc_qlen_PF_F_space(ledge1, best_t1, &ball_Fos);
    //IVP_DOUBLE qd1 = IVP_CLS.calc_qlen_PF_F_space(ledge1, best_t1, &ball_Fos);
    
    CORE;
    return IVP_FAULT;
}
#endif

#ifdef DEBUG
int situation(IVP_Cache_Ledge_Point *tri_cache, IVP_Cache_Ledge_Point *e_cache, IVP_Compact_Edge *triangle_e0, IVP_Compact_Edge *edge)
{
    // may be called from debugger (with print)
    
    // transforms triangle to flat (paper!)
    // also transforms edge
    // prints coords and labels

    // get world coords


    IVP_U_Point t_e_0, t_e_1, t_e_2;
    IVP_CLS.give_world_coords_AT(triangle_e0,tri_cache,&t_e_0);
    IVP_CLS.give_world_coords_AT(triangle_e0->get_next(),tri_cache,&t_e_1);
    IVP_CLS.give_world_coords_AT(triangle_e0->get_prev(),tri_cache,&t_e_2);


    IVP_U_Point e_0, e_1;
    IVP_CLS.give_world_coords_AT(edge,e_cache,&e_0);
    IVP_CLS.give_world_coords_AT(edge->get_next(),e_cache,&e_1);

    /*** search matrix flat on triangle with x-axis = triangle_e0 ***/

    // vector 0 (x) points in direction of triangle_e0
    IVP_U_Point vec_0;
    vec_0.subtract(&t_e_1, &t_e_0);
    vec_0.normize();

    // vector 2 (z) is orthogonal to triangle, i.e. emerging out of the paper
    IVP_U_Hesse vec_2; // world, normized
    IVP_CLS.calc_hesse_normized_AT( triangle_e0, tri_cache, &vec_2);
    vec_2.mult(-1.0f);
    
    // vector 1 (y) orthogonal to it
    IVP_U_Point vec_1;
    vec_1.calc_cross_product(&vec_0, &vec_2);
    vec_1.normize();

    // t_e_0 as origin
    IVP_U_Point trans;
    trans.set(&t_e_0);
    IVP_U_Matrix m_world_from_flat;
    m_world_from_flat.init_columns4(&vec_0, &vec_1, &vec_2, &trans);
    // m_flat_from_world.orthonormize();
    
    /*** transform whole situation ***/

    // transform triangle
    IVP_U_Point tf_0, tf_1, tf_2;
    m_world_from_flat.vimult4(&t_e_0, &tf_0);
    m_world_from_flat.vimult4(&t_e_1, &tf_1);
    m_world_from_flat.vimult4(&t_e_2, &tf_2);

    // transform edge
    IVP_U_Point ef_0, ef_1;
    m_world_from_flat.vimult4(&e_0, &ef_0);
    m_world_from_flat.vimult4(&e_1, &ef_1);

    // calc pierce position (if any)
    ivp_u_bool pierces = IVP_FALSE;
    ivp_u_bool pierce_inside = IVP_FALSE;
    IVP_U_Point pierce_point_f;

    IVP_U_Straight edge_straight;
    IVP_U_Point edge_vec;
    edge_vec.subtract(&e_1, &e_0);
    edge_straight.set(&e_0, &edge_vec);
    IVP_U_Hesse tri_hesse;
    IVP_CLS.calc_hesse_normized_AT( triangle_e0, tri_cache, &tri_hesse);

    IVP_U_Point pierce_point_w;
    IVP_DOUBLE q =0, r = 0;
    if( tri_hesse.calc_intersect_with(&edge_straight, &pierce_point_w) == IVP_OK ){
	pierces = IVP_TRUE;
	
	m_world_from_flat.vimult4(&pierce_point_w, &pierce_point_f);
	IVP_CLS.calc_qr_vals(triangle_e0,&pierce_point_w, &q, &r, tri_cache);
	if ( q >= 0 && r >= 0 && r+q<=1 ){
	    pierce_inside = IVP_TRUE;
	}
    }

    IVP_DOUBLE q_edge_start, q_edge_end, r_edge_start, r_edge_end;
    IVP_CLS.calc_qr_vals(triangle_e0, &e_0, &q_edge_start, &r_edge_start, tri_cache);
    IVP_CLS.calc_qr_vals(triangle_e0, &e_1, &q_edge_end, &r_edge_end, tri_cache);


    IVP_KK_Input kkin( triangle_e0, edge, tri_cache, e_cache);
    IVP_Unscaled_KK_Result kkr;
    IVP_CLS.calc_unscaled_KK_vals(kkin, &kkr);
    
    /** print everything ***/
    printf("triangle edge situation:");
    printf(" (coords so that triangle appears flat):\n\n");
    printf("Triangle: %f %f %f\n", tf_0.k[0], tf_0.k[1], tf_0.k[2]);
    printf("          %f %f %f\n", tf_1.k[0], tf_1.k[1], tf_1.k[2]);
    printf("          %f %f %f\n\n", tf_2.k[0], tf_2.k[1], tf_2.k[2]);
    printf("Edge:     %f %f %f\n", ef_0.k[0], ef_0.k[1], ef_0.k[2]);
    printf("          %f %f %f\n\n", ef_1.k[0], ef_1.k[1], ef_1.k[2]);
    
    if(pierces){
	printf("pierces %s triangle:  %f %f    q=%f r=%f \n",
	       (pierce_inside)?"inside":"outside",
	       pierce_point_f.k[0],
	       pierce_point_f.k[1],
	       q,
	       r);
    } else {
	printf("Edge is parallel.\n");
    }
    printf("\n");
    printf("edge startpoint q=%f r=%f\n", q_edge_start, r_edge_start);
    printf("                s=%f\n", IVP_CLS.calc_s_val(triangle_e0,&e_0, tri_cache));
    printf("     endpoint   q=%f r=%f\n", q_edge_end, r_edge_end);
    printf("                s=%f\n\n", IVP_CLS.calc_s_val(triangle_e0,&e_1, tri_cache));

    printf("tri_edge(K) to edge(L): slK=%f skL=%f\n",
	   kkr.checks_K[0] / (kkr.checks_K[0] + kkr.checks_K[1]) ,
	   kkr.checks_L[0] / (kkr.checks_L[0] + kkr.checks_L[1])
	   );
    return 0; // for debugger
}
#endif




///////////////////////////////////
///////////////////////////////////


#if 0

IVP_DOUBLE IVP_Mindist::calc_qlen_PF(IVP_Compact_Edge *P, IVP_Compact_Edge *F, IVP_Cache_Ledge_Point *m_cache_P, IVP_Cache_Ledge_Point *m_cache_F)
{
    IVP_U_Point *wP;
    wP = P->give_world_coords_AT(m_cache_P);

    IVP_DOUBLE q, r;
    F->calc_qr_vals(P, &q, &r, m_cache_F, m_cache_P);

    if ( q >= 0.0f && r >= 0.0f && r+q<=1.0f ){
	// inside triangle
	IVP_U_Hesse hesse_world;
	F->calc_hesse_AT(m_cache_F, &hesse_world);
	IVP_DOUBLE now_len = hesse_world.get_dist( wP );
	return now_len*now_len; // we want quad dists	
    }
    
    // outside triangle: take best edge instead
    IVP_DOUBLE checks[3];
    checks[0] = r;// go to direction if checks[j] < 0
    checks[1] = q;
    checks[2] = 1.0f - q - r;
    
    IVP_Compact_Edge *e;
    int j;
    IVP_DOUBLE minlen = P_DOUBLE_MAX;
    IVP_Compact_Edge *min_edge = 0;
    for (e=F,j=0;j<3;e=e->get_next(),j++){
	IVP_DOUBLE len = calc_qlen_PK( P, e, m_cache_P, m_cache_F);
	if(len < minlen){
	    min_edge = e;
	    minlen = len;
	}
    }
    IVP_ASSERT(min_edge); // what happened !?

    return minlen;
}




//////////////////////////////////
//////////////////////////////////

#ifdef DEBUG_CHECK_LEN
IVP_RETURN_TYPE IVP_Mindist::validate_termination_len(IVP_DOUBLE now_len)
{    
    IVP_DOUBLE t_len = IVP_Mindist::termination_len;

    // for debugger
    IVP_DOUBLE real_now_len = sqrt( now_len );
    IVP_DOUBLE real_term_len = sqrt( t_len );
    diff = real_now_len - real_term_len;
    if(now_len - IVP_MIN_TERMINATION_QDLEN > t_len){ // @@@ test for '=='
	CORE;
    }
    if(now_len - IVP_MIN_TERMINATION_QDLEN_EPS > t_len){ // @@@ test for '=='
	printf("Termination Len Violation: Epsilon Problem?\n");
	return IVP_FAULT;
    }
    if(now_len < ivp_mindist_settings.real_coll_dist * ivp_mindist_settings.real_coll_dist){
	CORE;
	return IVP_FAULT;
    }

    IVP_Mindist::termination_len = now_len;
    return IVP_OK;
}

#if 0
IVP_RETURN_TYPE IVP_Mindist::check_len_PP(IVP_Compact_Ledge *A, IVP_Compact_Edge *B, IVP_Cache_Ledge_Point *m_cache_A, IVP_Cache_Ledge_Point *m_cache_B)
{
    // termination function test
    IVP_DOUBLE now_len = calc_qlen_PP(A, B, m_cache_A, m_cache_B);
    return validate_termination_len(now_len);
}
#endif

IVP_RETURN_TYPE IVP_Mindist::check_len_PF(IVP_Compact_Edge *P, IVP_Compact_Edge *F, IVP_Cache_Ledge_Point *m_cache_P, IVP_Cache_Ledge_Point *m_cache_F)
{
    // termination function test
    IVP_DOUBLE now_len = calc_qlen_PF(P, F, m_cache_P, m_cache_F);

    return validate_termination_len(now_len);
}

IVP_RETURN_TYPE IVP_Mindist::check_len_PK(IVP_Compact_Edge *P, IVP_Compact_Edge *K, IVP_Cache_Ledge_Point *m_cache_P, IVP_Cache_Ledge_Point *m_cache_K)
{
    // termination function test
    IVP_DOUBLE now_len = calc_qlen_PK(P, K, m_cache_P, m_cache_K);
    return validate_termination_len(now_len);
}

IVP_RETURN_TYPE IVP_Mindist::check_len_KK(IVP_Compact_Edge *K, IVP_Compact_Edge *L, IVP_Cache_Ledge_Point *m_cache_K, IVP_Cache_Ledge_Point *m_cache_L)
{
    // termination function test
    IVP_DOUBLE now_len = calc_qlen_KK(K, L, m_cache_K, m_cache_L);

    return validate_termination_len(now_len);
}

IVP_RETURN_TYPE IVP_Mindist::check_len_FF(IVP_Compact_Edge *F, IVP_Compact_Edge *G, IVP_Cache_Ledge_Point *m_cache_F, IVP_Cache_Ledge_Point *m_cache_G)
{
    // termination function test
  
    IVP_DOUBLE now_len = IVP_CLS.p_optimize_FF(F, G, m_cache_F, m_cache_G, P_DOUBLE_MAX);

    return validate_termination_len(now_len);
}
#endif

#endif