|
◆ contng_simc()
integer function space_utilities::contng_simc |
( |
type (xy), dimension(:), intent(in) |
co, |
|
|
integer, intent(out) |
NT, |
|
|
integer, dimension(:), intent(out) |
IPT, |
|
|
integer, intent(out) |
NL, |
|
|
integer, dimension(:), intent(out) |
IPL |
|
) |
| |
|
private |
THIS SUBROUTINE PERFORMS TRIANGULATION.
IT DIVIDES THE X-Y PLANE INTO A NUMBER OF TRIANGLES ACCORDING TO GIVEN DATA POINTS IN THE PLANE, DETERMINES LINE SEGMENTS THAT FORM THE BORDER OF DATA AREA, AND DETERMINES THE TRIANGLE NUMBERS CORRESPONDING TO THE BORDER LINE SEGMENTS. AT COMPLETION, POINT NUMBERS OF THE VERTEXES OF EACH TRIANGLE ARE LISTED COUNTER-CLOCKWISE. POINT NUMBERS OF THE END POINTS OF EACH BORDER LINE SEGMENT ARE LISTED COUNTER-CLOCKWISE, LISTING ORDER OF THE LINE SEGMENTS BEING COUNTER-CLOCKWISE.
Return 0 if all right return 1 if IDENTICAL INPUT DATA POINTS FOUND return 2 if ALL DATA ARE COLLINEAR DATA POINTS - Parametri
-
[in] | co | ARRAY OF DIMENSION NDP CONTAINING THE COORDINATES OF THE DATA POINTS |
[out] | nt | NUMBER OF TRIANGLES |
[out] | nl | NUMBER OF BORDER LINE SEGMENTS |
[out] | ipt | ARRAY OF DIMENSION 6*NDP-15, WHERE THE POINT NUMBERS OF THE VERTEXES OF THE (IT)TH TRIANGLE ARE TO BE STORED AS THE (3*IT-2)ND, (3*IT-1)ST, AND (3*IT)TH ELEMENTS, IT=1,2,...,NT |
[out] | ipl | ARRAY OF DIMENSION 6*NDP, WHERE THE POINT NUMBERS OF THE END POINTS OF THE (IL)TH BORDER LINE SEGMENT AND ITS RESPECTIVE TRIANGLE NUMBER ARE TO BE STORED AS THE (3*IL-2)ND, (3*IL-1)ST, AND (3*IL)TH ELEMENTS, IL=1,2,..., NL. |
Definizione alla linea 371 del file space_utilities.F90.
380 call l4f_log(l4f_debug, "end triangulation first triangle")
396 dsqmn = dxmn**2+dymn**2
410 IF (ar.GE.(-armn) .AND. dsqi.GE.dsqmn) GO TO 230
417 230 ar = dy*dxmx-dx*dymx
418 IF (ar .LT. (-armx)) cycle
420 IF (ar.LE.armx .AND. dsqi.GE.dsqmx) cycle
427 IF (jpmx .LT. jpmn) jpmx = jpmx+nl0
437 ipl(jp3t3-2) = ipl(jp2t3-2)
438 ipl(jp3t3-1) = ipl(jp2t3-1)
439 ipl(jp3t3) = ipl(jp2t3)
443 ipl(jp2t3-2) = ipl(jp3t3-2)
444 ipl(jp2t3-1) = ipl(jp3t3-1)
445 ipl(jp2t3) = ipl(jp3t3)
456 l310 : DO jp2=jpmx,nl0
472 IF (jp2 == jpmx) then
480 ipl(nlnt3-1) = ipl(1)
489 IF (ipti.NE.ipl1 .AND. ipti.NE.ipl2) GO TO 300
491 IF (ipti.NE.ipl1 .AND. ipti.NE.ipl2) GO TO 300
496 300 IF (conxch_simc(co%X,co%Y,ip1,ipti,ipl1,ipl2) .EQ. 0) cycle l310
504 IF (jp2 .EQ. jpmx) ipl(jp2t3) = it
505 IF (jp2.EQ.nl0 .AND. ipl(3).EQ.it) ipl(3) = nt0
518 IF (nlf .EQ. 0) cycle l400
538 IF (ipl1.NE.ipt1 .AND. ipl1.NE.ipt2 .AND. ipl1.NE.ipt3) cycle
539 IF (ipl2.NE.ipt1 .AND. ipl2.NE.ipt2 .AND. ipl2.NE.ipt3) cycle
542 IF (ntf .EQ. 2) GO TO 330
544 IF (ntf .LT. 2) cycle l370
551 IF (ipti1.NE.ipl1 .AND. ipti1.NE.ipl2) GO TO 340
553 IF (ipti1.NE.ipl1 .AND. ipti1.NE.ipl2) GO TO 340
557 IF (ipti2.NE.ipl1 .AND. ipti2.NE.ipl2) GO TO 350
559 IF (ipti2.NE.ipl1 .AND. ipti2.NE.ipl2) GO TO 350
564 350 IF (conxch_simc(co%X,co%Y,ipti1,ipti2,ipl1,ipl2) .EQ. 0) cycle l370
589 IF ((iplj1.EQ.ipl1 .AND. iplj2.EQ.ipti2) .OR. &
590 (iplj2.EQ.ipl1 .AND. iplj1.EQ.ipti2)) ipl(jlt3) = itf(1)
591 IF ((iplj1.EQ.ipl2 .AND. iplj2.EQ.ipti1) .OR. &
592 (iplj2.EQ.ipl2 .AND. iplj1.EQ.ipti1)) ipl(jlt3) = itf(2)
597 IF (nlf .EQ. nlfc) cycle l400
604 DO jwl1=jwl1mn,nlft2,2
606 iwl(jwl-1) = iwl(jwl1-1)
613 call l4f_log(l4f_debug, "end triangulation appending")
623 IF (side(co(ip1),co(ip2),co(ip3)) .GE. 10.0) cycle
627 call l4f_log(l4f_debug, "end triangulation rearranging")
632 call l4f_log(l4f_debug, "end triangulation")
645 elemental double precision function vdsqf(co1,co2)
646 type(xy), intent(in) :: co1,co2
648 vdsqf = (co2%x-co1%x)**2+(co2%y-co1%y)**2
649 if (vdsqf == 0.d0) vdsqf = huge(vdsqf)
659 double precision function side(co1,co2,co3)
660 type(xy), intent(in):: co1,co2,co3
662 side = (co3%y-co1%y)*(co2%x-co1%x)-(co3%x-co1%x)*(co2%y-co1%y)
665 end function contng_simc
673 INTEGER FUNCTION conxch_simc (X,Y,I1,I2,I3,I4)
675 double precision, intent(in) :: X(:), Y(:)
678 integer, intent(in) :: I1,I2,I3,I4
680 double precision :: X0(4), Y0(4)
681 double precision :: C2SQ,C1SQ,A3SQ,B2SQ,B3SQ,A1SQ,A4SQ,B1SQ,B4SQ,A2SQ,C4SQ,C3SQ
683 double precision :: S1SQ,S2SQ,S3SQ,S4SQ,U1,U2,U3,U4
685 equivalence(c2sq,c1sq),(a3sq,b2sq),(b3sq,a1sq),(a4sq,b1sq), &
686 (b4sq,a2sq),(c4sq,c3sq)
697 u3 = (y0(2)-y0(3))*(x0(1)-x0(3))-(x0(2)-x0(3))*(y0(1)-y0(3))
698 u4 = (y0(1)-y0(4))*(x0(2)-x0(4))-(x0(1)-x0(4))*(y0(2)-y0(4))
699 IF (u3*u4 > 0.0) then
700 u1 = (y0(3)-y0(1))*(x0(4)-x0(1))-(x0(3)-x0(1))*(y0(4)-y0(1))
701 u2 = (y0(4)-y0(2))*(x0(3)-x0(2))-(x0(4)-x0(2))*(y0(3)-y0(2))
702 a1sq = (x0(1)-x0(3))**2+(y0(1)-y0(3))**2
703 b1sq = (x0(4)-x0(1))**2+(y0(4)-y0(1))**2
704 c1sq = (x0(3)-x0(4))**2+(y0(3)-y0(4))**2
705 a2sq = (x0(2)-x0(4))**2+(y0(2)-y0(4))**2
706 b2sq = (x0(3)-x0(2))**2+(y0(3)-y0(2))**2
707 c3sq = (x0(2)-x0(1))**2+(y0(2)-y0(1))**2
708 s1sq = u1*u1/(c1sq*max(a1sq,b1sq))
709 s2sq = u2*u2/(c2sq*max(a2sq,b2sq))
710 s3sq = u3*u3/(c3sq*max(a3sq,b3sq))
711 s4sq = u4*u4/(c4sq*max(a4sq,b4sq))
712 IF (min(s1sq,s2sq) < min(s3sq,s4sq)) idx = 1
717 END FUNCTION conxch_simc
Space utilities, derived from NCAR software.
|