[cig-commits] r16263 - seismo/2D/SPECFEM2D/trunk

[dkomati1 at geodynamics.org]

Author: dkomati1
Date: 2010-02-15 10:46:19 -0800 (Mon, 15 Feb 2010)
New Revision: 16263

Modified:
   seismo/2D/SPECFEM2D/trunk/specfem2D.F90
Log:
restored the sum for the rotation angle, added a comment explaining why,
and added a stop statement if anglerec is greater than zero and receivers are supposed to be oriented along the normal to the topography surface


Modified: seismo/2D/SPECFEM2D/trunk/specfem2D.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/specfem2D.F90	2010-02-15 12:36:51 UTC (rev 16262)
+++ seismo/2D/SPECFEM2D/trunk/specfem2D.F90	2010-02-15 18:46:19 UTC (rev 16263)
@@ -1859,7 +1859,6 @@
      allocate(adj_sourcearrays(1,1,1,1,1))
   endif
 
-
 if (ipass == 1) then
   if (nrecloc > 0) then
     allocate(anglerec_irec(nrecloc))
@@ -1872,6 +1871,10 @@
     allocate(sinrot_irec(1))
     allocate(rec_tangential_detection_curve(1))
   endif
+
+  if (rec_normal_to_surface .and. abs(anglerec) > 1.d-6) &
+    stop 'anglerec should be zero when receivers are normal to the topography'
+
   anglerec_irec(:) = anglerec * pi / 180.d0
   cosrot_irec(:) = cos(anglerec)
   sinrot_irec(:) = sin(anglerec)
@@ -1934,6 +1937,10 @@
 
     call tri_quad(n_tangential_detection_curve, n1_tangential_detection_curve, nnodes_tangential_curve)
 
+! in the case of a source force vector
+! users can give an angle with respect to the normal to the topography surface,
+! in which case we must compute the normal to the topography
+! and add it the existing rotation angle
     call compute_normal_vector( angleforce(i_source), nodes_tangential_curve(1,n_tangential_detection_curve(1)), &
       nodes_tangential_curve(1,n_tangential_detection_curve(2)), &
       nodes_tangential_curve(1,n_tangential_detection_curve(3)), &
@@ -7645,7 +7652,11 @@
       theta2 = - sign(1.d0,n3_x - n2_x) * acos(costheta2)
       theta3 = - sign(1.d0,n4_x - n3_x) * acos(costheta3)
 
-      angle = (theta1 + theta2 + theta3) / 3.d0 + PI/2.d0
+! a sum is needed here because in the case of a source force vector
+! users can give an angle with respect to the normal to the topography surface,
+! in which case we must compute the normal to the topography
+! and add it the existing rotation angle
+      angle = angle + (theta1 + theta2 + theta3) / 3.d0 + PI/2.d0
 
 end subroutine compute_normal_vector