Geometry/parallelepiped_point_dist_3d.f90 at master · lguez/Geometry · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
subroutine parallelepiped_point_dist_3d ( p1, p2, p3, p4, p, dist )

!*****************************************************************************80
!
!! PARALLELEPIPED_POINT_DIST_3D: distance ( parallelepiped, point ) in 3D.
!
!  Discussion:
!
!    A parallelepiped is a "slanted box", that is, opposite
!    sides are parallel planes.
!
!         *------------------*
!        / .                / \
!       /   .              /   \
!      /     .            /     \
!    P4------------------*       \
!      \        .         \       \
!       \        .         \       \
!        \        .         \       \
!         \       P2.........\.......\
!          \     .            \     /
!           \   .              \   /
!            \ .                \ /
!            P1-----------------P3
!
!  Licensing:
!
!    This code is distributed under the GNU LGPL license.
!
!  Modified:
!
!    03 January 2005
!
!  Author:
!
!    John Burkardt
!
!  Parameters:
!
!    Input, double precision P1(3), P2(3), P3(3), P4(3),
!    half of the corners of the box, from which the other corners can be
!    deduced.  The corners should be chosen so that the first corner
!    is directly connected to the other three.  The locations of
!    corners 5, 6, 7 and 8 will be computed by the parallelogram
!    relation.
!
!    Input, double precision P(3), the point which is to be checked.
!
!    Output, double precision DIST, the distance from the point to the box.
!    DIST is zero if the point lies exactly on the box.
!
  implicit none

  integer, parameter :: dim_num = 3

  double precision dis
  double precision dist
  double precision p(dim_num)
  double precision p1(dim_num)
  double precision p2(dim_num)
  double precision p3(dim_num)
  double precision p4(dim_num)
  double precision p5(dim_num)
  double precision p6(dim_num)
  double precision p7(dim_num)
  double precision p8(dim_num)
!
!  Fill in the other corners
!
  p5(1:dim_num) = p2(1:dim_num) + p3(1:dim_num) - p1(1:dim_num)
  p6(1:dim_num) = p2(1:dim_num) + p4(1:dim_num) - p1(1:dim_num)
  p7(1:dim_num) = p3(1:dim_num) + p4(1:dim_num) - p1(1:dim_num)
  p8(1:dim_num) = p2(1:dim_num) + p3(1:dim_num) + p4(1:dim_num) &
    - 2.0D+00 * p1(1:dim_num)
!
!  Compute the distance from the point P to each of the six
!  parallelogram faces.
!
  call parallelogram_point_dist_3d ( p1, p2, p3, p, dis )

  dist = dis

  call parallelogram_point_dist_3d ( p1, p2, p4, p, dis )

  dist = min ( dist, dis )

  call parallelogram_point_dist_3d ( p1, p3, p4, p, dis )

  dist = min ( dist, dis )

  call parallelogram_point_dist_3d ( p8, p5, p6, p, dis )

  dist = min ( dist, dis )

  call parallelogram_point_dist_3d ( p8, p5, p7, p, dis )

  dist = min ( dist, dis )

  call parallelogram_point_dist_3d ( p8, p6, p7, p, dis )

  dist = min ( dist, dis )

  return
end