[cig-commits] r16369 - in seismo/3D/SPECFEM3D_GLOBE/trunk: . UTILS/lib UTILS/seis_process
danielpeter at geodynamics.org
danielpeter at geodynamics.org
Tue Mar 2 14:29:51 PST 2010
Author: danielpeter
Date: 2010-03-02 14:29:50 -0800 (Tue, 02 Mar 2010)
New Revision: 16369
Added:
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/asc2sac.c
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/convolve_stf.c
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/phtimes.csh
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/readme.txt
Modified:
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/ascii2sac.csh
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_data.pl
seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_syn.pl
seismo/3D/SPECFEM3D_GLOBE/trunk/constants.h.in
seismo/3D/SPECFEM3D_GLOBE/trunk/meshfem3D.f90
seismo/3D/SPECFEM3D_GLOBE/trunk/specfem3D.f90
Log:
added a readme.txt to UTILS/seis_process and changed process_syn.pl & process_data.pl to call script phtimes.csh which uses a default binary provided by the iaspei-tau package (www.iris.edu) for phase picking
Added: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/asc2sac.c
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/asc2sac.c (rev 0)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/asc2sac.c 2010-03-02 22:29:50 UTC (rev 16369)
@@ -0,0 +1,131 @@
+/*
+ * asc2sac.c by Eh Tan
+ *
+ * Copyright (C) 2009, California Institute of Technology.
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+ *
+ */
+
+#include <stdio.h>
+#include <string.h>
+#include <stdlib.h>
+#include <limits.h>
+#include <errno.h>
+#include "sacio.h"
+#include "config.h"
+
+
+void asc2sac(char *ascfn) {
+ sac_int_t npts, max_npts, nerr, nconv, itmp;
+ float ftmp;
+ char *sacfn;
+ float *time, *data, a, b;
+ FILE *f;
+
+ sacfn = (char *)malloc(strlen(ascfn) + strlen(".sac") + 1);
+ if (!sacfn) {
+ fprintf(stderr, "Out of memory\n");
+ exit(1);
+ }
+ strcpy(sacfn, ascfn);
+ strcat(sacfn, ".sac");
+
+ /* reading ascii file */
+ f = fopen(ascfn, "r");
+ if(f == NULL) {
+ fprintf(stderr, "Cannot open file '%s' to read\n", ascfn);
+ exit(1);
+ }
+
+ time = data = 0;
+ max_npts = 0;
+ for (npts = 0; (nconv = fscanf(f, "%f %f\n", &a, &b)) == 2; ++npts) {
+ if (nconv != 2) {
+ fprintf(stderr, "error while reading file '%s'\n", ascfn);
+ exit(1);
+ }
+ if (npts >= max_npts) {
+ max_npts = max_npts ? 2 * max_npts : 1024;
+ time = (float*) realloc(time, max_npts * sizeof(float));
+ data = (float*) realloc(data, max_npts * sizeof(float));
+ if(time == NULL || data == NULL) {
+ fprintf(stderr, "Out of memory\n");
+ exit(1);
+ }
+ }
+ time[npts] = a;
+ data[npts] = b;
+ }
+ if (nconv != EOF || ferror(f)) {
+ fprintf(stderr, "error while reading file '%s' (on or near line %ld)\n", ascfn, (long)npts);
+ exit(1);
+ }
+ fclose(f);
+ /* finished reading ascii file */
+
+ /* write SAC data usng SAC IO library */
+ nerr = 0;
+ newhdr();
+ setnhv("npts", &npts, &nerr, strlen("npts"));
+ itmp = 1;
+ setlhv("leven", &itmp, &nerr, strlen("leven"));
+ ftmp = time[1] - time[0];
+ setfhv("delta", &ftmp, &nerr, strlen("delta"));
+ setfhv("b", &(time[0]), &nerr, strlen("b"));
+ setfhv("e", &(time[npts-1]), &nerr, strlen("e"));
+ setihv("iftype", "itime", &nerr, strlen("iftype"), strlen("itime"));
+ setihv("idep", "idisp", &nerr, strlen("idep"), strlen("idisp"));
+
+ if(nerr) {
+ fprintf(stderr, "error when setting header for '%s'\n", sacfn);
+ exit(1);
+ }
+
+ wsac0(sacfn, time, data, &nerr, strlen(sacfn));
+
+ if(nerr) {
+ fprintf(stderr, "error when writing '%s'\n", sacfn);
+ exit(1);
+ }
+
+ free(time);
+ free(data);
+ free(sacfn);
+
+ return;
+}
+
+
+int main(int argc, char *argv[]) {
+ int i;
+
+ if (argc < 2) {
+ fprintf(stderr,
+ "usage: %s FILE...\n"
+ "\n"
+ " Converts ASCII time series files to SAC binary files.\n"
+ " Input files must have two columns: t, f(t).\n"
+ "\n",
+ argv[0]);
+ exit(1);
+ }
+
+ for (i = 1; i < argc; ++i) {
+ asc2sac(argv[i]);
+ }
+
+ return 0;
+}
Added: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/convolve_stf.c
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/convolve_stf.c (rev 0)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/lib/convolve_stf.c 2010-03-02 22:29:50 UTC (rev 16369)
@@ -0,0 +1,270 @@
+/*
+ * convolve_stf.c by Qinya Liu, Eh Tan
+ *
+ * Copyright (C) 2003-2009, California Institute of Technology.
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or
+ * (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+ *
+ */
+
+#include <stdio.h>
+#include <string.h>
+#include <stdlib.h>
+#include <errno.h>
+#include <math.h>
+#include "sacio.h"
+#include "config.h"
+
+
+/* defined in libsac */
+void fft(float *xreal, float *ximag, sac_int_t n, sac_int_t idir);
+
+
+static void fzero(float *dst, sac_int_t n) {
+ while (n)
+ dst[--n] = 0.0;
+}
+
+static void fcpy(float *dst, const float *src, sac_int_t n) {
+ while (n) {
+ --n;
+ dst[n] = src[n];
+ }
+}
+
+void convolve(float **pconv, sac_int_t *pnconv,
+ const float *data, sac_int_t ndata,
+ const float *stf, sac_int_t nstf)
+{
+ sac_int_t nconv, ncorr, i;
+ struct { float *xreal, *ximag; } cdata, cstf, ccorr;
+ float *conv, *buffer;
+
+ nconv = ndata + nstf - 1;
+ conv = (float *)malloc(nconv * sizeof(float));
+
+ for (ncorr = 2; ncorr < nconv; ncorr *= 2)
+ ;
+
+ buffer = (float *)malloc(2 * 3 * ncorr * sizeof(float));
+ cdata.xreal = buffer + 0 * ncorr;
+ cdata.ximag = buffer + 1 * ncorr;
+ cstf.xreal = buffer + 2 * ncorr;
+ cstf.ximag = buffer + 3 * ncorr;
+ ccorr.xreal = buffer + 4 * ncorr;
+ ccorr.ximag = buffer + 5 * ncorr;
+
+ fcpy(cdata.xreal, data, ndata);
+ fzero(cdata.xreal + ndata, ncorr - ndata);
+ fzero(cdata.ximag, ncorr);
+
+ fcpy(cstf.xreal, stf, nstf);
+ fzero(cstf.xreal + nstf, ncorr - nstf);
+ fzero(cstf.ximag, ncorr);
+
+ fft(cdata.xreal, cdata.ximag, ncorr, 1);
+ fft(cstf.xreal, cstf.ximag, ncorr, 1);
+
+ for (i = 0; i < ncorr; ++i) {
+ /* ccorr[i] = cdata[i] * cstf[i] */
+ ccorr.xreal[i] = cdata.xreal[i] * cstf.xreal[i] - cdata.ximag[i] * cstf.ximag[i];
+ ccorr.ximag[i] = cdata.xreal[i] * cstf.ximag[i] + cdata.ximag[i] * cstf.xreal[i];
+ }
+
+ fft(ccorr.xreal, ccorr.ximag, ncorr, -1);
+ for (i = 0; i < nconv; ++i) {
+ conv[i] = ccorr.xreal[i] * ncorr;
+ }
+
+ free(buffer);
+
+ *pconv = conv;
+ *pnconv = nconv;
+
+ return;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+ char cstf, *endpt;
+ float hdur, *data;
+ int j;
+ sac_int_t datasize;
+ const int min_nhdur = 10;
+
+ if(argc < 4) {
+ fprintf(stderr,
+ "Usage: %s g[|t] HDUR FILE...\n"
+ " This program convolves SAC files with a gaussian(|triangle)\n"
+ " source time function of given half duration\n",
+ argv[0]);
+ return 1;
+ }
+
+ if(strcmp(argv[1], "t") != 0 && strcmp(argv[1], "g") != 0) {
+ fprintf(stderr,"The source time function type could only be triangle(t) or gaussian(g) \n");
+ return 1;
+ }
+ cstf = argv[1][0];
+
+ errno = 0;
+ hdur = (float)strtod(argv[2], &endpt);
+ if(errno || endpt == argv[2] || *endpt != 0) {
+ fprintf(stderr,"No floating point number can be formed from %s\n",argv[2]);
+ return 1;
+ }
+
+
+ /* loop over input sac files */
+ datasize = 0;
+ data = NULL;
+ for (j=3; j<argc; j++) {
+ sac_int_t max, npts, nlen, nerr;
+ float beg, del, dt, origin, tmp[1];
+ sac_int_t nstf, nconv, i;
+ float hstf, *stf, *conv;
+ char *outf;
+
+
+ fprintf(stderr, "convolving sac file %s with half duration %7.3f\n",
+ argv[j], hdur);
+
+ if((outf = (char *) malloc(strlen(argv[j]) + strlen(".conv") + 1)) == NULL) {
+ fprintf(stderr,"Out of memory\n");
+ return 1;
+ }
+ strcpy(outf, argv[j]);
+ strcat(outf, ".conv");
+
+ /* read header to get the length of time series */
+ max = 1;
+ rsac1(argv[j], tmp, &nlen, &beg, &del, &max, &nerr, strlen(argv[j]));
+ if(nerr != -803) {
+ fprintf(stderr,"Error reading sac file %s\n",argv[j]);
+ return 1;
+ }
+
+ getnhv("npts", &npts, &nerr, strlen("npts"));
+
+ /* now we know how much memory we need to allocate */
+ if(npts > datasize) {
+ data = (float*)realloc(data, npts * sizeof(float));
+ if(data == NULL) {
+ fprintf(stderr, "out of memory\n");
+ return 1;
+ }
+ datasize = npts;
+ }
+
+ /* read the data */
+ max = npts;
+ rsac1(argv[j], data, &nlen, &beg, &del, &max, &nerr, strlen(argv[j]));
+ if(nerr) {
+ fprintf(stderr,"Error reading sac file %s\n",argv[j]);
+ return 1;
+ }
+
+ /* get additional info */
+ getfhv("delta", &dt, &nerr, strlen("delta"));
+
+ getfhv("o", &origin, &nerr, strlen("o"));
+ if(nerr) {
+ /* Assuming origin time is 0, per convention of SPECFEM */
+ fprintf(stderr, " Ignoring the error, assuming origin time is 0\n");
+ origin = 0.0;
+ setfhv("o", &origin, &nerr, strlen("o"));
+ }
+
+
+
+ /* creat source time function time series */
+ if(min_nhdur * dt / 2 > hdur) {
+ fprintf(stderr,"The half duration %f is too small to convolve\n", hdur);
+ return 1;
+ }
+
+ nstf = (int)ceil(2 * hdur/dt)+1;
+ hstf = (nstf-1)*dt/2;
+
+ if((stf = (float *) malloc(nstf*sizeof(float))) == NULL) {
+ fprintf(stderr,"Error in allocating memory for source time function\n");
+ return 1;
+ }
+ if(cstf == 't') {
+ /* triangular */
+ for (i=0; i<nstf/2; i++)
+ stf[i] = i * dt / (hstf * hstf);
+ for (i=nstf/2; i<nstf; i++)
+ stf[i] = (2 * hstf - i * dt)/ (hstf*hstf);
+ } else {
+ /* gaussian */
+ const float decay_rate = 1.628;
+ float alpha = decay_rate / hdur;
+ float divisor = sqrt(4*atan(1.0) /*pi*/);
+ for (i=0; i<nstf; i++) {
+ float tao_i = fabs(i*dt - hstf);
+ stf[i] = alpha * exp(- alpha*alpha* tao_i*tao_i) / divisor;
+ }
+ }
+
+
+ /* creat convolution time series */
+ convolve(&conv,&nconv,data,npts,stf,nstf);
+ for(i=0; i<nconv; i++)
+ conv[i] *= dt;
+
+
+ /* update sac file header */
+ {
+ int n;
+ float bnew, enew, aminm, amaxm, amean;
+
+ bnew = origin - hstf + beg;
+ enew = bnew + del * nconv;
+
+ aminm = amaxm = conv[0];
+ amean = 0;
+ for(n=0; n<nconv; n++) {
+ aminm = (aminm < conv[n]) ? aminm : conv[n];
+ amaxm = (amaxm > conv[n]) ? amaxm : conv[n];
+ amean += conv[n];
+ }
+ amean /= nconv;
+
+ setfhv("b", &bnew, &nerr, strlen("b"));
+ setfhv("e", &enew, &nerr, strlen("e"));
+ setnhv("npts", &nconv, &nerr, strlen("npts"));
+ setfhv("depmin", &aminm, &nerr, strlen("depmin"));
+ setfhv("depmax", &amaxm, &nerr, strlen("depmax"));
+ setfhv("depmen", &amean, &nerr, strlen("depmen"));
+ }
+
+ /* output to .conv sac file */
+ wsac0(outf, data, conv, &nerr, strlen(outf));
+
+ if(nerr) {
+ fprintf(stderr, "Not able to write to file %s\n", argv[j]);
+ return 1;
+ }
+
+ free(conv);
+ free(stf);
+ free(outf);
+ }
+
+ return 0;
+}
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/ascii2sac.csh
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/ascii2sac.csh 2010-03-02 20:54:18 UTC (rev 16368)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/ascii2sac.csh 2010-03-02 22:29:50 UTC (rev 16369)
@@ -1,6 +1,9 @@
#!/bin/csh -f
-
+#
# Vala Hjorleifsdottir and Qinya Liu, Caltech, Jan 2007
+#
+# uses `asc2sac` binary which can be compiled from UTILS/lib/asc2sac.c
+# (requires SAC libraries for compilation)
foreach file ($*)
echo $file
Added: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/phtimes.csh
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/phtimes.csh (rev 0)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/phtimes.csh 2010-03-02 22:29:50 UTC (rev 16369)
@@ -0,0 +1,60 @@
+#!/bin/csh
+#
+# script to output travel times
+# based on iaspei-tau package
+# intended to work for P and S phases
+#
+# usage: ./phtimes.csh depth(km) distance(degrees) phase
+#
+# e.g. ./phtimes.csh 611 52.474 P
+# ./phtimes.csh 611 52.474 S
+#
+###########################################
+# MODEL PARAMETER
+
+# locates your IASPEI-TAU installation
+set dir="/opt/seismo-util/source/iaspei-tau"
+
+# uses IASP91 tables
+set model="$dir/tables/iasp91"
+
+###########################################
+
+# input arguments
+set depth="$1"
+set delta="$2"
+set branch="$3"
+
+# checks arguments
+if( $depth =~ "" ) then
+ echo "Usage: phtimes.csh depth(km) distance(degrees) phase"
+ exit 1
+endif
+if( $delta =~ "" ) then
+ echo "Usage: phtimes.csh depth(km) distance(degrees) phase"
+ exit 1
+endif
+if( $branch =~ "" ) then
+ echo "Usage: phtimes.csh depth(km) distance(degrees) phase"
+ exit 1
+endif
+
+
+# travel times
+$dir/bin/ttimes <<EOF > tmp.log
+Y
+$model
+$branch
+
+N
+$depth
+$delta
+-1
+-1
+EOF
+
+# outputs single line with phase name and traveltime
+fgrep -A 1 "#" tmp.log | tail -n 1 | awk '{print $3,$4}'
+
+
+# in case of problems, see the created tmp.log for further informations
\ No newline at end of file
Property changes on: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/phtimes.csh
___________________________________________________________________
Name: svn:executable
+ *
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_data.pl
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_data.pl 2010-03-02 20:54:18 UTC (rev 16368)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_data.pl 2010-03-02 22:29:50 UTC (rev 16369)
@@ -76,7 +76,7 @@
if ($opt_f and not ($opt_R or $opt_i)) {die("-f option goes together with -i or -R\n");}
$saclst="saclst";
-$phtimes="phtimes";
+$phtimes="./phtimes.csh";
$eps=1e-5; $undef=-12345.0; $cundef="-12345";
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_syn.pl
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_syn.pl 2010-03-02 20:54:18 UTC (rev 16368)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/process_syn.pl 2010-03-02 22:29:50 UTC (rev 16369)
@@ -70,8 +70,9 @@
$eps = 0.1;
$saclst="saclst";
-$phtimes="phtimes";
+$phtimes="./phtimes.csh";
$asc2sac="ascii2sac.csh";
+
#if (! -e $saclst) {die(" No $saclst file\n");}
#if (! -e $phtimes) {die("No $phtimes file\n");}
#if (! -e $asc2sac) {die("No $asc2sac file\n");}
Added: seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/readme.txt
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/readme.txt (rev 0)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/UTILS/seis_process/readme.txt 2010-03-02 22:29:50 UTC (rev 16369)
@@ -0,0 +1,99 @@
+--------------------------------
+readme
+--------------------------------
+
+To compare synthetics and data, the following steps are recommended:
+
+1. Make sure that both synthetic and observed seismograms have the
+ correct station/event and timing information.
+
+2. Convolve synthetic seismograms with a source time function with the
+ half duration specified in the CMTSOLUTION file, provided, as recommended,
+ you used a zero half duration in the SEM simulations.
+
+3. Resample both observed and synthetic seismograms to a common sampling rate.
+
+4. Cut the records using the same window.
+
+5. Remove the trend and mean from the records and taper them.
+
+6. Remove the instrument response from the observed seismograms (recommended)
+ or convolve the synthetic seismograms with the instrument response.
+
+7. Make sure that you apply the same filters to both observed and synthetic seismograms.
+ Preferably, avoid filtering your records more than once.
+
+8. Now, you are ready to compare your synthetic and observed seismograms.
+
+
+Example:
+--------
+
+- data processing:
+
+ > process_data.pl -m CMTSOLUTION -s 1.0 -l 0/4000 -i -f -t 40/500 -p -x bp DATA/1999.330*.LH?.SAC
+
+ - adds CMTSOLUTION information
+ - resamples ( 1 s sampling rate)
+ - cuts ( between 0 and 4000 s)
+ - transfers data to displacement and filters (between 40 and 500 s)
+ - picks P & S arrivals
+ - adds suffix 'bp'
+
+- synthetics processing:
+
+ > process_syn.pl -m CMTSOLUTION -h -a STATIONS -s 1.0 -l 0/4000 -f -t 40/500 -p -x bp SEM/*.LH?.sem
+
+ - adds CMTSOLUTION information
+ - convolves with source-time function
+ - adds STATIONS information
+ - resamples ( 1 s sampling rate)
+ - cuts ( between 0 and 4000 s)
+ - filters (between 40 and 500 s) compatible with transfer function
+ - picks P & S arrivals
+ - adds suffix 'bp'
+
+to rotate the horizontal components:
+
+ > rotate.pl -l 0 -L 4000 -d DATA/*.LHE.SAC.bp
+ > rotate.pl -l 0 -L 4000 SEM/.LHE.semd.sac.bp
+
+
+
+Note:
+-----
+
+the processing scripts make use of additional software
+(assumed to be installed in a default directory /opt/seismo-util/):
+
+- SAC
+ http://www.iris.washington.edu/software/sac/
+
+ dowload software package (IRIS DMC Software - Processing Programs)
+
+ the scripts 'process_**.pl' make use of the provided binary `sac` and `saclst`.
+
+- iaspei-tau
+ http://www.iris.washington.edu/pub/programs/iaspei-tau/
+
+ download software package (IRIS DMC Software - Processing Programs)
+ and install it by running:
+ > cd iaspei-tau/
+ > ./build-package
+
+ the scripts 'process_**.pl' call 'phtimes.csh' (for picking phases)
+ which makes use of the provided binary `ttimes`.
+ please modify the paths in the script 'phtimes.csh' accordingly.
+
+- asc2sac
+ source code provided in UTILS/lib or on https://geodynamics.org/portals/seismo/ inside
+ the post-process utilities:
+ https://geodynamics.org/portals/seismo/samples/seis_process-1.0.1.tar.gz
+
+- convolve_stf
+ source code provided in UTILS/lib or on https://geodynamics.org/portals/seismo/ inside
+ the post-process utilities:
+ https://geodynamics.org/portals/seismo/samples/seis_process-1.0.1.tar.gz
+
+
+
\ No newline at end of file
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/constants.h.in
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/constants.h.in 2010-03-02 20:54:18 UTC (rev 16368)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/constants.h.in 2010-03-02 22:29:50 UTC (rev 16369)
@@ -110,7 +110,7 @@
!! pathname of the topography file (un-smoothed)
! character (len=*), parameter :: PATHNAME_TOPO_FILE = 'DATA/topo_bathy/ETOPO1.xyz'
-! Use GLL points to capture TOPOGRAPHY and ELLIPTICITY
+! Use GLL points to capture TOPOGRAPHY and ELLIPTICITY (experimental feature)
logical,parameter :: USE_GLL = .false.
! maximum depth of the oceans in trenches and height of topo in mountains
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/meshfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/meshfem3D.f90 2010-03-02 20:54:18 UTC (rev 16368)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/meshfem3D.f90 2010-03-02 22:29:50 UTC (rev 16369)
@@ -410,6 +410,10 @@
! be mainly a performance reason. changing the codes to adopt modules
! will have to prove that it performs as fast as it does without now.
!
+! another reason why modules are avoided, is to make the code thread safe.
+! having different threads access the same data structure and modifying it at the same time
+! would lead to problems. passing arguments is a way to avoid such complications.
+!
! however, the mesher makes one exception here: it uses the
! module "meshfem3D_models_par" defined in the 'meshfem3D_models.f90' file.
! the exception is based on the fact, that when one wants to incorporate
@@ -440,7 +444,7 @@
!
! finally, in order to compile the new mesher with your new file(s),
! you will add it to the list in the 'Makefile.in' file and run
-! ./configure to recreate a new Makefile.
+! `configure` to recreate a new Makefile.
!
!
!-------------------------------------------------------------------------------------------------
Modified: seismo/3D/SPECFEM3D_GLOBE/trunk/specfem3D.f90
===================================================================
--- seismo/3D/SPECFEM3D_GLOBE/trunk/specfem3D.f90 2010-03-02 20:54:18 UTC (rev 16368)
+++ seismo/3D/SPECFEM3D_GLOBE/trunk/specfem3D.f90 2010-03-02 22:29:50 UTC (rev 16369)
@@ -788,7 +788,12 @@
! passing them along as arguments to the routine makes the code slower.
! it seems that this stack/heap criterion is more complicated.
!
-! also, it would nice to test, if using modules with static arrays would perform as well as this.
+! another reason why modules are avoided, is to make the code thread safe.
+! having different threads access the same data structure and modifying it at the same time
+! would lead to problems. passing arguments is a way to avoid such complications.
+!
+! nevertheless, it would be nice to test - where possible - , if using modules
+! together with static arrays would perform as well as this.
! at least, it would make the code more elegant and less error prone...
!
! note 2: in general, most of the computation time for our earthquake simulations is spent
More information about the CIG-COMMITS
mailing list