[CIG-SEISMO] Noise simulation in SPECFEM2D

Ryan Modrak rmodrak at gmail.com
Mon Jan 28 20:14:06 PST 2013


Hi Ying,

Thank you for your interest in the noise tomography paper and code.

Here are some comments to add to what Yang said.

If you haven't already, download SPECFEM2D, run the configuration script,
and run EXAMPLES/noise_uniform/process.sh.  I just did so using the latest
release (SPECFEM2D-7.0.0.tar.gz<http://www.geodynamics.org/cig/software/specfem2d/SPECFEM2D-7.0.0.tar.gz>).
When I plot the resulting density kernels (contained in CAB_POS/proc* and
CBA_NEG/proc*), I get

     [image: Inline image 2]     [image: Inline image 1]


The above figures match the ones in Tromp et al. 2010, but if you look very
closely you can see subtle lines that run opposite to the main sensitivity
contours.  These lines are artefacts produced by imperfect absorbing
boundaries. We could eliminate these artefacts by enlarging the model while
keeping the simulation time the same.

A similar sort of effect occurs with the cross correlations.  Because the
boundaries are imperfect, a part of the generating wavefield is reflected
back from the boundary, which together with the placement of the master
receiver closer to one side than the other, causes the two cross
correlation branches to have different amplitudes.  If you or I were to
extend the boundaries far enough, there would be no boundary interactions,
and the two branches would have the same amplitude.

These issues are trivial to Yang and me because we've spent many hours
writing and debugging the code, but given the lack of detailed discussion
in the paper (it being long enough already as it is), I would not expect
the underlying rationale for, say, the sign change you pointed out, to be
in any way obvious. As we publish more about noise tomography and implement
new features (such as PML), I can only hope the clarity of the presentation
improves.

Please don't hesitate to email me if you have any other questions.

Best regards,
Ryan



On Mon, Jan 28, 2013 at 9:34 PM, Ying Ying <yingzi.ying at me.com> wrote:

> Hi Yang,
>
> Thank you so much.
>
>  I am always happy to eliminate TRIVIAL bugs to make things come to
> obvious. Your suggestions are really helpful.  Sincerely obliged to your
> warm help.
>
> Best,
> Yingzi
>
>
> Sent from my iPad
>
> On Jan 29, 2013, at 2:12, Yang Luo <yangl.princeton at gmail.com> wrote:
>
> Your result seems good to me (the asymmetry of the two branches may be due
> to the fact that the computational domain is slightly smaller than the one
> used for the paper).
> Note that the time shown on your plot is what I referred to as the
> "simulation time", which is always positive (0 to 120 s), whereas the time
> shown in the paper is the "cross-correlation time", which spans -60 to 60 s
> (in this case, "simulation time"="cross-correlation time"+60 ).
> You may see that for the "cross-correlatioin time", pulse in the positive
> branch
> appears around 33 s, which corresponds to the 100 km distance (wavespeed
> is 3 km/s)
>
> the minus sign only flips the waveform, which is trivial.
>
> the feature in 3D (double the NSTEP) is for real observations.
> Imagine you have noisy seismograms of 100 s at each station,
> then the cross correlation between each pair of them leads to 200 s cross
> correlations. The feature makes it easier for user to put the length of
> noisy seismograms directly, instead of the length of the cross correlations.
> It seems to me that feature is NOT included in 2D code,
> but it is trivial too.
>
>
>
>
>
>
> On Mon, Jan 28, 2013 at 3:56 PM, Yingzi Ying <yingzi.ying at me.com> wrote:
>
>> Hi Yang,
>>
>> I checked the code "noise_tomography.f90" again and attached a snip bellow
>> ------------------------------**-------------------------
>>   elseif( time_function_type == 4 ) then
>>     !reproduce time function from Figure 2a of Tromp et al. 2010
>>     do it = 1,NSTEP
>>       t = it*deltat
>>       time_function_noise(it) = factor_noise * &
>>        4.*aval**2. * (3. - 12.*aval*(t-t0)**2. + 4.*aval**2.*(t-t0)**4.)
>> * &
>>        exp(-aval*(t-t0)**2.)
>>     enddo
>>
>> ------------------------------**-------------------------
>>
>> Is there should be a "-" sign before "factor_noise"? If I did this, the
>> result will be just corrected, even the length of cross-correlation is
>> still in controversial.
>>
>>
>> Best,
>> Yingzi
>> On 01/28/2013 08:54 PM, Yang Luo wrote:
>>
>>>  Hi Jeroen & Yingzi,
>>> I didn't implement the 2D version, which was actually done by Ryan
>>> (cc-ed) based on the 3D code we have already had.
>>> But I just checked the 2D code again and ran the example Yingzi
>>> mentioned, all looks good to me.
>>> The code I have may not be current, since I cannot use svn to get the
>>> latest version.
>>> But since the older version I have works, I assume the latest one should
>>> work as well.
>>> Yingzi, please see my following two cents, and make sure that you didn't
>>> misunderstand what you have already obtained from your tests.
>>> 1. the positive branch --- it refers to the second half (i.e., 60 to 120
>>> s for the example) of the seismograms you obtained (not for t>0 in the
>>> simulation).
>>> In other words, our simulation actually should start from t=-60s, but to
>>> avoid introducing this complexity, the time axis has been shifted.
>>> Given this into consideration, I think what you saw (i.e., two envelops)
>>> in your seismograms is correct.
>>> You just missed the point that the time from the simulation (almost
>>> always positive) is not the cross-correlation time (could be negative).
>>> Tromp et. al. (2010) provides step-by-step figures, which I think you
>>> might not fully understand yet.
>>> 2. the time_function_type --- you are right, the default should be
>>> changed to the one that reads in the file S_squares.
>>> But for benchmarking purposes, it is fine to put the one we have used in
>>> our paper.
>>> (Ryan, could you please fix that? just provide a default S_squres file
>>> there using the time_function you hard coded).
>>> 3. nu_master --- it refers to the noise direction, which is trivial in
>>> 2D.
>>> Yingzi, should you still be confused or claim that what you saw is
>>> totally different from the paper,
>>> please provide a picture of what you saw.
>>> Thanks,
>>> Yang
>>>
>>>
>>> On Mon, Jan 28, 2013 at 11:45 AM, Jeroen Tromp <jtromp at princeton.edu
>>>  <mailto:jtromp at princeton.edu>> wrote:
>>>
>>>     Hi Yang:
>>>
>>>     A question for you.
>>>
>>>     Best regards,
>>>
>>>     Jeroen
>>>
>>>
>>>     Begin forwarded message:
>>>
>>>      *From:* Yingzi Ying <yingzi.ying at me.com <mailto:yingzi.ying at me.com
>>>> >>
>>>>     *Date:* January 28, 2013, 11:18:00 PST
>>>>     *To:* "cig-seismo at geodynamics.org
>>>>     <mailto:cig-seismo@**geodynamics.org <cig-seismo at geodynamics.org>>"
>>>> <cig-seismo at geodynamics.org
>>>>     <mailto:cig-seismo@**geodynamics.org <cig-seismo at geodynamics.org>>>
>>>>     *Subject:* *[CIG-SEISMO] Noise simulation in SPECFEM2D*
>>>>     *Reply-To:* <cig-seismo at geodynamics.org
>>>>     <mailto:cig-seismo@**geodynamics.org <cig-seismo at geodynamics.org>>>
>>>>
>>>>
>>>>
>>>>     Dear Tromp, Dimitri and all,
>>>>
>>>>     I am interested in noise simulation with SPECFEM2D(V7.0.0) to get
>>>>     ensemble cross-correlation. I tried to reproduce the result in
>>>>     fig. 2(a)
>>>>     in Tromp's paper "Noise cross-correlation sensitivity kernels" by
>>>>     executing ./process.sh in EXAMPLES/noise_uniform/.
>>>>
>>>>     My calculated interferogram after step 2 has two envelopes in the
>>>>     positive branch, and their positions are unconcerned with the time
>>>>     delay
>>>>     between two stations. This is totally different with the Tromp's
>>>>     result.
>>>>
>>>>     I checked the source code "noise_tomography.f90" and found the value
>>>>     "time_function_type" is fixed to "4", which means the variable
>>>>     "time_function_type" in "DATA/SOURCE" will not take action and the
>>>>     "DATA/NOISE_TOMOGRAPHY/S_**squared" will not be read. Furthermore,
>>>>     there
>>>>     is no "nu_master" in SPECFEM2D, which exists in 3D version.
>>>>
>>>>     I don't know what goes wrong, and what I should do to get a correct
>>>>     ensemble interferogram. Anyone who knows please kindly give me a
>>>> hint.
>>>>
>>>>     Thanks in advance.
>>>>
>>>>     Best,
>>>>     Yingzi
>>>>     ______________________________**_________________
>>>>     CIG-SEISMO mailing list
>>>>     CIG-SEISMO at geodynamics.org <mailto:CIG-SEISMO@**geodynamics.org<CIG-SEISMO at geodynamics.org>
>>>> >
>>>>     http://geodynamics.org/cgi-**bin/mailman/listinfo/cig-**seismo<http://geodynamics.org/cgi-bin/mailman/listinfo/cig-seismo>
>>>>
>>>
>>>
>>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://geodynamics.org/pipermail/cig-seismo/attachments/20130128/3fd7d3b3/attachment-0001.htm 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: image/png
Size: 57741 bytes
Desc: not available
Url : http://geodynamics.org/pipermail/cig-seismo/attachments/20130128/3fd7d3b3/attachment-0002.png 
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: image/png
Size: 39919 bytes
Desc: not available
Url : http://geodynamics.org/pipermail/cig-seismo/attachments/20130128/3fd7d3b3/attachment-0003.png 


More information about the CIG-SEISMO mailing list