diff options
Diffstat (limited to 'beam_tracing/octave/beam2face_distance.m')
| -rw-r--r-- | beam_tracing/octave/beam2face_distance.m | 47 |
1 files changed, 47 insertions, 0 deletions
diff --git a/beam_tracing/octave/beam2face_distance.m b/beam_tracing/octave/beam2face_distance.m new file mode 100644 index 0000000..d7d5a61 --- /dev/null +++ b/beam_tracing/octave/beam2face_distance.m @@ -0,0 +1,47 @@ +function [hit_distance, hit_position, is_face_hit] = beam2face_distance(beam,face) + %% return distance to face if beam hits it or infinity otherwise + % for beams and faces structure description see face_beam_interaction.m + + % default returning values for the case when beam misses the face + is_face_hit = false; + hit_position = [NA, NA]; + hit_distance = Inf; + + k = beam.k; + % face direction or kf + kf=face.vertex2 - face.vertex1; % not a unit vector + + %% simple check for intersection of two vectors + % if beam are parallel no intersection is possible + % we do this via simulated cross product calculation + angle_tolerance = 1e-15; + if ( abs( (k(1)*kf(2)-k(2)*kf(1))/norm(k)/norm(kf) ) <= angle_tolerance ) + % beams never intercept + return; + end + + %% let's find the intersection of lines passing through face and light beam + % we introduce parameters tb and tf + % which define distance (in arb. units) along k vectors + % we are solving + % xb_o+ k_x*tb = v1_x+ kf_x*tf + % yb_o+ k_y*tb = v1_y+ kf_y*tf + % A*t = B + A = [ k(1) , -kf(1); k(2), -kf(2)]; + B= [ face.vertex1(1)-beam.origin(1); face.vertex1(2)-beam.origin(2) ]; + t = A\B; tb=t(1); tf=t(2); + + % beam intersects face only if 0<=tf<=1 and tb>=0 + if ( (0 <= tf) && (tf<=1) && tb >=0 ) + % intersection, beam hits the face + is_face_hit = true; + else + % beam misses the face + return; + end + + % calculating hit position + hit_position=face.vertex1+kf*tf; + hit_distance = tb; % distance along the light beam +end + |