Gaussian source #279
Gaussian source #279
Conversation
There was a problem hiding this comment.
Awesome, with such a source we can build now Mehdis CDM source non-analytically, but it would be aware of velocity model etc. Which opens us the door for the Volcanology guys. Of cours we would need to talk with Mehdi about if it still would be valid...
Please find my comments below...
| nucleation_x=None, nucleation_y=None, | ||
| tref=0.0, decimation_factor=1): | ||
| def two_d_gaussian(x, y, amplitude, xo, yo, sigma_x, sigma_y, theta_g, offset): | ||
| xo = float(xo) |
There was a problem hiding this comment.
It provides much more functionality than you are using for now. Could you please add a docstring to the input of this function? What is theta_g for. For the purpose of future extension.
There was a problem hiding this comment.
Cool stuff! theta_g should turn the Gaussian moment/slip ellipse away from the x and y axis. Nice feature. Maybe one can bifurcate to the more simple and cheap expression if theta_g=0.
However, if the ellipse is turned, then the tail ends of the distribution may stick out of the plane
(same for x0 ne 0 or y0 ne 0). Would this compromise the moment of the source? We put in the amplitude=1 but may get out less.
| @@ -200,32 +221,51 @@ def discretize_rect_source(deltas, deltat, north, east, depth, | |||
| dist = num.sqrt(dist_x**2 + dist_y**2) | |||
| times = dist / velocity | |||
|
|
|||
| if sigma_x is not None and sigma_y is not None: | |||
| slip_function = two_d_gaussian( | |||
There was a problem hiding this comment.
Now you leave the Gaussian always centered. Maybe tie the offset to the nuc_x, nuc_y attributes or better introduce new ones as the slip-maximum doesnt need to be at the hypocentre?! But sure that adds the 73rd argument ....
Here you set theta always to zero. I dont understand what theta is. But we also have to thing in terms of performance- if you leave it always zero the equation reduces to:
(1 / (2 * num.pi * sigma ** 2)) * num.exp( - ((gridpoints - gp0) ** 2) / (2 * sigma ** 2))
And we could save the expensive calculation of the sin- and cos-terms? Simply because this function is called millions of times and it might add up-what do you think? Maybe simply an if statement in case both are zero?
There was a problem hiding this comment.
I think the slip distribution would in parts leave the rupture plane if y0, x0 ne 0. As for a rotated slip distribution, I am unsure if this messes up the moment of the source.
| @@ -1692,20 +1740,38 @@ class RectangularSource(SourceWithDerivedMagnitude): | |||
| optional=True, | |||
| default=1) | |||
|
|
|||
| gaussian_cutoff = Float.T( | |||
| optional=True, | |||
| help='if given, turn on 2D gaussian distributed slip on source ' | |||
There was a problem hiding this comment.
From this description I still dont get what the cut-off here does. I would need more explanation to know how to specify it.
| @@ -1692,20 +1740,38 @@ class RectangularSource(SourceWithDerivedMagnitude): | |||
| optional=True, | |||
| default=1) | |||
|
|
|||
There was a problem hiding this comment.
As mentioned above maybe already finalize the attributes and allow the centre of the gaussian to be optionally shifting as well...The rest of the code allows it already. Now you simply hardcode it to zero.
RectangularSource, activated if source parametergaussian_cutoffis setdepth_minanddepth_max)ruputure_plane_rake, angle to rotate the orientation of the rectangle in the source plane.