[cig-commits] r18133 - seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts
carltape at geodynamics.org
carltape at geodynamics.org
Thu Mar 24 11:36:29 PDT 2011
Author: carltape
Date: 2011-03-24 11:36:29 -0700 (Thu, 24 Mar 2011)
New Revision: 18133
Added:
seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xyz_to_vtk.pl
Removed:
seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xy_to_vtk.pl
Modified:
seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/split_sr_vtk.pl
Log:
rename convert_xy_to_vtk.pl --> convert_xyz_to_vtk.pl and added comments; modified split_sr_vtk.pl to allow finite source model
Deleted: seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xy_to_vtk.pl
===================================================================
--- seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xy_to_vtk.pl 2011-03-24 15:46:29 UTC (rev 18132)
+++ seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xy_to_vtk.pl 2011-03-24 18:36:29 UTC (rev 18133)
@@ -1,51 +0,0 @@
-#!/usr/bin/perl -w
-
-# this program converts the xy file (format as in psxy, can be segmented
-# by >) to the vtk polydata format to be viewed in Paraview.
-# Qinya Liu, May 2007, Caltech
-
-if (@ARGV == 0) { die("Usage: convert_xy_to_vtk.pl xyfile\n");}
-
-$xyfile=$ARGV[0];
-$vtkfile="${xyfile}.vtk";
-
-open(FILE,">$vtkfile");
-
-print FILE "# vtk DataFile Version 2.0\n";
-print FILE "Really cool data\n";
-print FILE "ASCII\n";
-print FILE "DATASET POLYDATA\n";
-
-# figure out the total number of points
-open(XY,"$xyfile") || die("Check $xyfile\n");
- at xy = <XY>;
-close(XY);
-# record the points
-$total_points = 0; $points="";
-# record the lines
-$nlines = 0; $npoints = 0; $nnum = 0; $lines=""; $lines2="";
-
-for ($i=0;$i<@xy;$i++) {
- @tmp = split(" ",$xy[$i]);
- if ($tmp[0] ne ">") {
- if (@tmp == 2) {
- $points.="$tmp[0] $tmp[1] 0.0\n";}
- elsif (@tmp == 3) {
- $points.="$tmp[0] $tmp[1] $tmp[2]\n";}
- else {die("Check this line @tmp\n");}
- $lines.="$total_points ";
- $total_points++; $npoints++;
- }
- if ($tmp[0] eq ">" or $i == @xy-1) {
- # change to a new line
- $nlines ++ ;
- $lines2.="$npoints $lines\n";
- $nnum += $npoints+1;
- $npoints=0; $lines="";
- }
-}
-print FILE "POINTS $total_points float\n";
-print FILE "$points";
-print FILE "LINES $nlines $nnum\n";
-print FILE "$lines2";
-close(FILE);
Copied: seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xyz_to_vtk.pl (from rev 18132, seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xy_to_vtk.pl)
===================================================================
--- seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xyz_to_vtk.pl (rev 0)
+++ seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/convert_xyz_to_vtk.pl 2011-03-24 18:36:29 UTC (rev 18133)
@@ -0,0 +1,53 @@
+#!/usr/bin/perl -w
+
+# convert_xyz_to_vtk.pl
+# This program converts an xy file or xyz file to a vtk file.
+# It will also work with segmented files with ">" as the delimiter.
+# Qinya Liu, May 2007, Caltech
+
+if (@ARGV == 0) { die("Usage: convert_xyz_to_vtk.pl xyzfile\n");}
+
+$xyzfile=$ARGV[0];
+$vtkfile="${xyzfile}.vtk";
+
+open(FILE,">$vtkfile");
+
+print FILE "# vtk DataFile Version 2.0\n";
+print FILE "XYZ data\n";
+print FILE "ASCII\n";
+print FILE "DATASET POLYDATA\n";
+
+# figure out the total number of points
+open(XYZ,"$xyzfile") || die("Check $xyzfile\n");
+ at xyz = <XYZ>;
+close(XYZ);
+# record the points
+$total_points = 0; $points="";
+# record the lines
+$nlines = 0; $npoints = 0; $nnum = 0; $lines=""; $lines2="";
+
+for ($i=0;$i<@xyz;$i++) {
+ @tmp = split(" ",$xyz[$i]);
+ if ($tmp[0] ne ">") {
+ # if input has 2 columns, set 3rd column to zero (flat earth scenario, z up, e.g. UTM)
+ if (@tmp == 2) {
+ $points.="$tmp[0] $tmp[1] 0.0\n";}
+ elsif (@tmp == 3) {
+ $points.="$tmp[0] $tmp[1] $tmp[2]\n";}
+ else {die("Check this line @tmp\n");}
+ $lines.="$total_points ";
+ $total_points++; $npoints++;
+ }
+ if ($tmp[0] eq ">" or $i == @xyz-1) {
+ # change to a new line
+ $nlines ++ ;
+ $lines2.="$npoints $lines\n";
+ $nnum += $npoints+1;
+ $npoints=0; $lines="";
+ }
+}
+print FILE "POINTS $total_points float\n";
+print FILE "$points";
+print FILE "LINES $nlines $nnum\n";
+print FILE "$lines2";
+close(FILE);
Modified: seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/split_sr_vtk.pl
===================================================================
--- seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/split_sr_vtk.pl 2011-03-24 15:46:29 UTC (rev 18132)
+++ seismo/3D/SPECFEM3D/trunk/utils/Visualization/vtk_scripts/split_sr_vtk.pl 2011-03-24 18:36:29 UTC (rev 18133)
@@ -8,19 +8,25 @@
# epicenter.vtk
# receiver.vtk
#
-# Future modifications:
-# 1: allow for multiple sources (currently assumes ONE only)
-# 2: compute epicenter for spherical model (assumes FLAT model only)
+# future modifications:
+# compute epicenter for spherical/elliptical model (assumes FLAT model only),
+# or simply have SPECFEM3D output the three files instead of sr.vtk
#
# EXAMPLE:
-# split_sr_vtk.pl sr.vtk .
-# split_sr_vtk.pl sr.vtk
+# split_sr_vtk.pl sr_tmp.vtk . 1
#
#-------------------
-if (@ARGV == 0) { die("Usage: split_sr_vtk.pl srfile odir\n");}
+if (@ARGV < 3) { die("Usage: split_sr_vtk.pl srfile odir nsrc\n");}
$srfile = $ARGV[0];
+$odir = $ARGV[1];
+$nsrc = $ARGV[2];
+
+print "src-rec file is $srfile\n";
+print "output directory for files is $odir\n";
+print "number of sources is $nsrc (see CMTSOLUTION)\n";
+
if($ARGV[1]) {$odir = $ARGV[1]} else {$odir = "."}
#$odir = $ARGV[0];
if(not -e $odir) {$odir = "."}
@@ -35,38 +41,40 @@
($xsrc,$ysrc,undef) = split(" ",$lines[5]);
$nlines = @lines;
-$nrec = $nlines - 7;
+$nrec = $nlines - $nsrc - 6;
-# write source file (hypocenter)
-open(FILE1,">$odir/source.vtk");
+print "nlines = $nlines, nsrc = $nsrc , nrec = $nrec\n";
+
+# write receiver file
+open(FILE1,">$odir/receiver.vtk");
print FILE1 "$lines[0]";
-print FILE1 "Source VTK file\n";
+print FILE1 "Receiver VTK file\n";
print FILE1 "$lines[2]";
print FILE1 "$lines[3]";
-print FILE1 "POINTS 1 float\n";
-print FILE1 "$lines[5]";
+print FILE1 "POINTS $nrec float\n";
+print FILE1 "@lines[$nsrc+5..$nrec+$nsrc+4]";
close(FILE1);
-# write source file (epicenter)
-$zepi = "0.0";
-open(FILE2,">$odir/epicenter.vtk");
+# write source file (hypocenter)
+open(FILE2,">$odir/source.vtk");
print FILE2 "$lines[0]";
-print FILE2 "Epicenter VTK file\n";
+print FILE2 "Source VTK file\n";
print FILE2 "$lines[2]";
print FILE2 "$lines[3]";
-print FILE2 "POINTS 1 float\n";
-print FILE2 "$xsrc $ysrc $zepi";
+print FILE2 "POINTS $nsrc float\n";
+print FILE2 "@lines[5..$nsrc+4]";
close(FILE2);
-# write receiver file
+# write source file (epicenter)
+# note: needs to be revised for non-planar surfaces (e.g., SPECFEM3D_GLOBE)
$zepi = "0.0";
-open(FILE3,">$odir/receiver.vtk");
+open(FILE3,">$odir/epicenter.vtk");
print FILE3 "$lines[0]";
-print FILE3 "Receiver VTK file\n";
+print FILE3 "Epicenter VTK file\n";
print FILE3 "$lines[2]";
print FILE3 "$lines[3]";
-print FILE3 "POINTS $nrec float\n";
-print FILE3 "@lines[6..$nlines-1]";
+print FILE3 "POINTS 1 float\n";
+print FILE3 "$xsrc $ysrc $zepi";
close(FILE3);
print "\n-- writing out source.vtk, epicenter.vtk, receiver.vtk in directory $odir --\n";
More information about the CIG-COMMITS
mailing list