C++ cross_v3_v3v3函数代码示例

static void get_strand_normal(Material *ma, const float surfnor[3], float surfdist, float nor[3])
{
	float cross[3], nstrand[3], vnor[3], blend;

	if (!((ma->mode & MA_STR_SURFDIFF) || (ma->strand_surfnor > 0.0f)))
		return;

	if (ma->mode & MA_STR_SURFDIFF) {
		cross_v3_v3v3(cross, surfnor, nor);
		cross_v3_v3v3(nstrand, nor, cross);

		blend = dot_v3v3(nstrand, surfnor);
		CLAMP(blend, 0.0f, 1.0f);

		interp_v3_v3v3(vnor, nstrand, surfnor, blend);
		normalize_v3(vnor);
	}
	else {
		copy_v3_v3(vnor, nor);
	}

	if (ma->strand_surfnor > 0.0f) {
		if (ma->strand_surfnor > surfdist) {
			blend = (ma->strand_surfnor - surfdist) / ma->strand_surfnor;
			interp_v3_v3v3(vnor, vnor, surfnor, blend);
			normalize_v3(vnor);
		}
	}

	copy_v3_v3(nor, vnor);
}

/* our own triangle intersection, so we can fully control the epsilons and
 * prevent corner case from going wrong*/
static int meshdeform_tri_intersect(const float orig[3], const float end[3], const float vert0[3],
                                    const float vert1[3], const float vert2[3],
                                    float r_isectco[3], float r_uvw[3])
{
	float edge1[3], edge2[3], tvec[3], pvec[3], qvec[3];
	float det, inv_det, u, v, dir[3], isectdir[3];

	sub_v3_v3v3(dir, end, orig);

	/* find vectors for two edges sharing vert0 */
	sub_v3_v3v3(edge1, vert1, vert0);
	sub_v3_v3v3(edge2, vert2, vert0);

	/* begin calculating determinant - also used to calculate U parameter */
	cross_v3_v3v3(pvec, dir, edge2);

	/* if determinant is near zero, ray lies in plane of triangle */
	det = dot_v3v3(edge1, pvec);

	if (UNLIKELY(det == 0.0f)) {
		return 0;
	}

	inv_det = 1.0f / det;

	/* calculate distance from vert0 to ray origin */
	sub_v3_v3v3(tvec, orig, vert0);

	/* calculate U parameter and test bounds */
	u = dot_v3v3(tvec, pvec) * inv_det;
	if (u < -EPSILON || u > 1.0f + EPSILON)
		return 0;

	/* prepare to test V parameter */
	cross_v3_v3v3(qvec, tvec, edge1);

	/* calculate V parameter and test bounds */
	v = dot_v3v3(dir, qvec) * inv_det;
	if (v < -EPSILON || u + v > 1.0f + EPSILON)
		return 0;

	r_isectco[0] = (1.0f - u - v) * vert0[0] + u * vert1[0] + v * vert2[0];
	r_isectco[1] = (1.0f - u - v) * vert0[1] + u * vert1[1] + v * vert2[1];
	r_isectco[2] = (1.0f - u - v) * vert0[2] + u * vert1[2] + v * vert2[2];

	r_uvw[0] = 1.0f - u - v;
	r_uvw[1] = u;
	r_uvw[2] = v;

	/* check if it is within the length of the line segment */
	sub_v3_v3v3(isectdir, r_isectco, orig);

	if (dot_v3v3(dir, isectdir) < -EPSILON)
		return 0;
	
	if (dot_v3v3(dir, dir) + EPSILON < dot_v3v3(isectdir, isectdir))
		return 0;

	return 1;
}

/**
 * Add randomness to stroke
 * \param gps: Stroke data
 * \param brush: Brush data
 */
void gp_randomize_stroke(bGPDstroke *gps, bGPDbrush *brush)
{
	bGPDspoint *pt1, *pt2, *pt3;
	float v1[3];
	float v2[3];
	if (gps->totpoints < 3) {
		return;
	}

	/* get two vectors using 3 points */
	pt1 = &gps->points[0];
	pt2 = &gps->points[1];
	pt3 = &gps->points[(int)(gps->totpoints * 0.75)];

	sub_v3_v3v3(v1, &pt2->x, &pt1->x);
	sub_v3_v3v3(v2, &pt3->x, &pt2->x);
	normalize_v3(v1);
	normalize_v3(v2);

	/* get normal vector to plane created by two vectors */
	float normal[3];
	cross_v3_v3v3(normal, v1, v2);
	normalize_v3(normal);

	/* get orthogonal vector to plane to rotate random effect */
	float ortho[3];
	cross_v3_v3v3(ortho, v1, normal);
	normalize_v3(ortho);

	/* Read all points and apply shift vector (first and last point not modified) */
	for (int i = 1; i < gps->totpoints - 1; ++i) {
		bGPDspoint *pt = &gps->points[i];
		/* get vector with shift (apply a division because random is too sensitive */
		const float fac = BLI_frand() * (brush->draw_random_sub / 10.0f);
		float svec[3];
		copy_v3_v3(svec, ortho);
		if (BLI_frand() > 0.5f) {
			mul_v3_fl(svec, -fac);
		}
		else {
			mul_v3_fl(svec, fac);
		}

		/* apply shift */
		add_v3_v3(&pt->x, svec);
	}

}

float angle_signed_on_axis_v3v3v3_v3(const float v1[3], const float v2[3], const float v3[3], const float axis[3])
{
	float v1_proj[3], v2_proj[3], tproj[3];
	float angle;

	sub_v3_v3v3(v1_proj, v1, v2);
	sub_v3_v3v3(v2_proj, v3, v2);

	/* project the vectors onto the axis */
	project_v3_v3v3(tproj, v1_proj, axis);
	sub_v3_v3(v1_proj, tproj);

	project_v3_v3v3(tproj, v2_proj, axis);
	sub_v3_v3(v2_proj, tproj);

	angle = angle_v3v3(v1_proj, v2_proj);

	/* calculate the sign (reuse 'tproj') */
	cross_v3_v3v3(tproj, v2_proj, v1_proj);
	if (dot_v3v3(tproj, axis) < 0.0f) {
		angle = ((float)(M_PI * 2.0)) - angle;
	}

	return angle;
}

/* adjust bone roll to align Z axis with vector
 * vec is in local space and is normalized
 */
float ED_rollBoneToVector(EditBone *bone, const float align_axis[3], const short axis_only)
{
	float mat[3][3], nor[3];

	sub_v3_v3v3(nor, bone->tail, bone->head);
	vec_roll_to_mat3(nor, 0.0f, mat);
	
	/* check the bone isn't aligned with the axis */
	if (!is_zero_v3(align_axis) && angle_v3v3(align_axis, mat[2]) > FLT_EPSILON) {
		float vec[3], align_axis_proj[3], roll;
		
		/* project the new_up_axis along the normal */
		project_v3_v3v3(vec, align_axis, nor);
		sub_v3_v3v3(align_axis_proj, align_axis, vec);
		
		if (axis_only) {
			if (angle_v3v3(align_axis_proj, mat[2]) > (float)(M_PI / 2.0)) {
				negate_v3(align_axis_proj);
			}
		}
		
		roll = angle_v3v3(align_axis_proj, mat[2]);
		
		cross_v3_v3v3(vec, mat[2], align_axis_proj);
		
		if (dot_v3v3(vec, nor) < 0) {
			roll = -roll;
		}
		
		return roll;
	}

	return 0.0f;
}

/**
 * finalize after accumulation.
 */
static void calc_tangent_ortho(float ts[3][3])
{
	float v_tan_a[3], v_tan_b[3];
	float t_vec_a[3], t_vec_b[3];

	normalize_v3(ts[2]);

	copy_v3_v3(v_tan_a, ts[0]);
	copy_v3_v3(v_tan_b, ts[1]);

	cross_v3_v3v3(ts[1], ts[2], v_tan_a);
	mul_v3_fl(ts[1], dot_v3v3(ts[1], v_tan_b) < 0.0f ? -1.0f : 1.0f);

	/* orthognalise tangent */
	mul_v3_v3fl(t_vec_a, ts[2], dot_v3v3(ts[2], v_tan_a));
	sub_v3_v3v3(ts[0], v_tan_a, t_vec_a);

	/* orthognalise bitangent */
	mul_v3_v3fl(t_vec_a, ts[2], dot_v3v3(ts[2], ts[1]));
	mul_v3_v3fl(t_vec_b, ts[0], dot_v3v3(ts[0], ts[1]) / dot_v3v3(v_tan_a, v_tan_a));
	sub_v3_v3(ts[1], t_vec_a);
	sub_v3_v3(ts[1], t_vec_b);

	normalize_v3(ts[0]);
	normalize_v3(ts[1]);
}

/* adjust bone roll to align Z axis with vector
 * vec is in local space and is normalized
 */
float ED_armature_ebone_roll_to_vector(const EditBone *bone, const float align_axis[3], const bool axis_only)
{
	float mat[3][3], nor[3];
	float vec[3], align_axis_proj[3], roll = 0.0f;

	BLI_ASSERT_UNIT_V3(align_axis);

	sub_v3_v3v3(nor, bone->tail, bone->head);

	/* If tail == head or the bone is aligned with the axis... */
	if (normalize_v3(nor) <= FLT_EPSILON || (fabsf(dot_v3v3(align_axis, nor)) >= (1.0f - FLT_EPSILON))) {
		return roll;
	}

	vec_roll_to_mat3_normalized(nor, 0.0f, mat);

	/* project the new_up_axis along the normal */
	project_v3_v3v3_normalized(vec, align_axis, nor);
	sub_v3_v3v3(align_axis_proj, align_axis, vec);

	if (axis_only) {
		if (angle_v3v3(align_axis_proj, mat[2]) > (float)(M_PI_2)) {
			negate_v3(align_axis_proj);
		}
	}

	roll = angle_v3v3(align_axis_proj, mat[2]);

	cross_v3_v3v3(vec, mat[2], align_axis_proj);

	if (dot_v3v3(vec, nor) < 0.0f) {
		return -roll;
	}
	return roll;
}

/**
 * Return true if the 2x axis are both aligned when projected into the view.
 * In this case, we can't usefully project the cursor onto the plane.
 */
static bool isPlaneProjectionViewAligned(const TransInfo *t)
{
  const float eps = 0.001f;
  const float *constraint_vector[2];
  int n = 0;
  for (int i = 0; i < 3; i++) {
    if (t->con.mode & (CON_AXIS0 << i)) {
      constraint_vector[n++] = t->con.mtx[i];
      if (n == 2) {
        break;
      }
    }
  }
  BLI_assert(n == 2);

  float view_to_plane[3], plane_normal[3];

  getViewVector(t, t->center_global, view_to_plane);

  cross_v3_v3v3(plane_normal, constraint_vector[0], constraint_vector[1]);
  normalize_v3(plane_normal);

  float factor = dot_v3v3(plane_normal, view_to_plane);
  return fabsf(factor) < eps;
}

static void testAxialSymmetry(BGraph *graph, BNode *root_node, BNode *node1, BNode *node2, BArc *arc1, BArc *arc2, float axis[3], float limit, int group)
{
	const float limit_sq = limit * limit;
	float nor[3], vec[3], p[3];

	sub_v3_v3v3(p, node1->p, root_node->p);
	cross_v3_v3v3(nor, p, axis);

	sub_v3_v3v3(p, root_node->p, node2->p);
	cross_v3_v3v3(vec, p, axis);
	add_v3_v3(vec, nor);
	
	cross_v3_v3v3(nor, vec, axis);
	
	if (fabsf(nor[0]) > fabsf(nor[1]) && fabsf(nor[0]) > fabsf(nor[2]) && nor[0] < 0) {
		negate_v3(nor);
	}
	else if (fabsf(nor[1]) > fabsf(nor[0]) && fabsf(nor[1]) > fabsf(nor[2]) && nor[1] < 0) {
		negate_v3(nor);
	}
	else if (fabsf(nor[2]) > fabsf(nor[1]) && fabsf(nor[2]) > fabsf(nor[0]) && nor[2] < 0) {
		negate_v3(nor);
	}
	
	/* mirror node2 along axis */
	copy_v3_v3(p, node2->p);
	BLI_mirrorAlongAxis(p, root_node->p, nor);
	
	/* check if it's within limit before continuing */
	if (len_squared_v3v3(node1->p, p) <= limit_sq) {
		/* mark node as symmetric physically */
		copy_v3_v3(root_node->symmetry_axis, nor);
		root_node->symmetry_flag |= SYM_PHYSICAL;
		root_node->symmetry_flag |= SYM_AXIAL;

		/* flag side on arcs */
		flagAxialSymmetry(root_node, node1, arc1, group);
		flagAxialSymmetry(root_node, node2, arc2, group);
		
		if (graph->axial_symmetry) {
			graph->axial_symmetry(root_node, node1, node2, arc1, arc2);
		}
	}
	else {
		/* NOT SYMMETRIC */
	}
}

static bool convex(const float p0[3], const float up[3], const float a[3], const float b[3])
{
	/* Vec3 va = a-p0, vb = b-p0; */
	float va[3], vb[3], tmp[3];
	sub_v3_v3v3(va, a, p0);
	sub_v3_v3v3(vb, b, p0);
	cross_v3_v3v3(tmp, va, vb);
	return dot_v3v3(up, tmp) >= 0;
}

static int convex(float *p0, float *up, float *a, float *b)
{
	// Vec3 va = a-p0, vb = b-p0;
	float va[3], vb[3], tmp[3];
	sub_v3_v3v3(va, a, p0);
	sub_v3_v3v3(vb, b, p0);
	cross_v3_v3v3(tmp, va, vb);
	return dot_v3v3(up, tmp) >= 0;
}

static bool depth_read_normal(
        const ViewContext *vc, const bglMats *mats, const int mval[2],
        float r_normal[3])
{
	/* pixels surrounding */
	bool  depths_valid[9] = {false};
	float coords[9][3] = {{0}};

	ARegion *ar = vc->ar;
	const ViewDepths *depths = vc->rv3d->depths;

	for (int x = 0, i = 0; x < 2; x++) {
		for (int y = 0; y < 2; y++) {
			const int mval_ofs[2] = {mval[0] + (x - 1), mval[1] + (y - 1)};

			const double depth = (double)depth_read_zbuf(vc, mval_ofs[0], mval_ofs[1]);
			if ((depth > depths->depth_range[0]) && (depth < depths->depth_range[1])) {
				if (depth_unproject(ar, mats, mval_ofs, depth, coords[i])) {
					depths_valid[i] = true;
				}
			}
			i++;
		}
	}

	const int edges[2][6][2] = {
	    /* x edges */
	    {{0, 1}, {1, 2},
	     {3, 4}, {4, 5},
	     {6, 7}, {7, 8}},
	    /* y edges */
	    {{0, 3}, {3, 6},
	     {1, 4}, {4, 7},
	     {2, 5}, {5, 8}},
	};

	float cross[2][3] = {{0.0f}};

	for (int i = 0; i < 6; i++) {
		for (int axis = 0; axis < 2; axis++) {
			if (depths_valid[edges[axis][i][0]] && depths_valid[edges[axis][i][1]]) {
				float delta[3];
				sub_v3_v3v3(delta, coords[edges[axis][i][0]], coords[edges[axis][i][1]]);
				add_v3_v3(cross[axis], delta);
			}
		}
	}

	cross_v3_v3v3(r_normal, cross[0], cross[1]);

	if (normalize_v3(r_normal) != 0.0f) {
		return true;
	}
	else {
		return false;
	}
}

static void flush_pixel(const MResolvePixelData *data, const int x, const int y)
{
	float st[2] = {(x + 0.5f) / data->w, (y + 0.5f) / data->h};
	float *st0, *st1, *st2;
	float *tang0, *tang1, *tang2;
	float no0[3], no1[3], no2[3];
	float fUV[2], from_tang[3][3], to_tang[3][3];
	float u, v, w, sign;
	int r;

	const int i0 = data->i0;
	const int i1 = data->i1;
	const int i2 = data->i2;

	st0 = data->mtface[data->face_index].uv[i0];
	st1 = data->mtface[data->face_index].uv[i1];
	st2 = data->mtface[data->face_index].uv[i2];

	multiresbake_get_normal(data, no0, data->face_index, i0);   /* can optimize these 3 into one call */
	multiresbake_get_normal(data, no1, data->face_index, i1);
	multiresbake_get_normal(data, no2, data->face_index, i2);

	resolve_tri_uv(fUV, st, st0, st1, st2);

	u = fUV[0];
	v = fUV[1];
	w = 1 - u - v;

	if (data->pvtangent) {
		tang0 = data->pvtangent + data->face_index * 16 + i0 * 4;
		tang1 = data->pvtangent + data->face_index * 16 + i1 * 4;
		tang2 = data->pvtangent + data->face_index * 16 + i2 * 4;

		/* the sign is the same at all face vertices for any non degenerate face.
		 * Just in case we clamp the interpolated value though. */
		sign = (tang0[3] * u + tang1[3] * v + tang2[3] * w) < 0 ? (-1.0f) : 1.0f;

		/* this sequence of math is designed specifically as is with great care
		 * to be compatible with our shader. Please don't change without good reason. */
		for (r = 0; r < 3; r++) {
			from_tang[0][r] = tang0[r] * u + tang1[r] * v + tang2[r] * w;
			from_tang[2][r] = no0[r] * u + no1[r] * v + no2[r] * w;
		}

		cross_v3_v3v3(from_tang[1], from_tang[2], from_tang[0]);  /* B = sign * cross(N, T)  */
		mul_v3_fl(from_tang[1], sign);
		invert_m3_m3(to_tang, from_tang);
	}
	else {
		zero_m3(to_tang);
	}

	data->pass_data(data->lores_dm, data->hires_dm, data->bake_data,
	                data->ibuf, data->face_index, data->lvl, st, to_tang, x, y);
}

/**
 * accumulate edge-vectors from all polys.
 */
static void calc_tangent_loop_accum(const float v_dir_prev[3],
                                    const float v_dir_next[3],
                                    float r_tspace[3][3])
{
  add_v3_v3v3(r_tspace[1], v_dir_prev, v_dir_next);

  if (compare_v3v3(v_dir_prev, v_dir_next, FLT_EPSILON * 10.0f) == false) {
    const float weight = fabsf(acosf(dot_v3v3(v_dir_next, v_dir_prev)));
    float nor[3];

    cross_v3_v3v3(nor, v_dir_prev, v_dir_next);
    normalize_v3(nor);

    cross_v3_v3v3(r_tspace[0], r_tspace[1], nor);

    mul_v3_fl(nor, weight);
    /* accumulate weighted normals */
    add_v3_v3(r_tspace[2], nor);
  }
}

/* builds an X and a Y axis from the given Z axis */
static void build_coordinate_frame(float axisX[3], float axisY[3], const float axisZ[3])
{
	const float faX = fabsf(axisZ[0]);
	const float faY = fabsf(axisZ[1]);
	const float faZ = fabsf(axisZ[2]);

	if (faX <= faY && faX <= faZ) {
		const float len = sqrtf(axisZ[1] * axisZ[1] + axisZ[2] * axisZ[2]);
		axisY[0] = 0; axisY[1] = axisZ[2] / len; axisY[2] = -axisZ[1] / len;
		cross_v3_v3v3(axisX, axisY, axisZ);
	}
	else if (faY <= faZ) {
		const float len = sqrtf(axisZ[0] * axisZ[0] + axisZ[2] * axisZ[2]);
		axisX[0] = axisZ[2] / len; axisX[1] = 0; axisX[2] = -axisZ[0] / len;
		cross_v3_v3v3(axisY, axisZ, axisX);
	}
	else {
		const float len = sqrtf(axisZ[0] * axisZ[0] + axisZ[1] * axisZ[1]);
		axisX[0] = axisZ[1] / len; axisX[1] = -axisZ[0] / len; axisX[2] = 0;
		cross_v3_v3v3(axisY, axisZ, axisX);
	}
}