How can I extract data with geographic coordinate system from vtk file?

ASPECT
1

Hello,

I have a problem how to extract geographic coordinates from the solution pvtk file.

When I tried to perform the simulation using 3D box geometry, I could get cartesian coordinates using paraview, and I could convert them to geographic coordinates (because I converted input parameters from geographic coordinates to cartesian before the simulations, I could convert to geographic coordinates again).

However, I tried to simulate using 3D chunk, It needed r (radius), phi (longitude), theta (90-latitude). I expected that the coordinates of the simulation results were also r, phi, theta coordinates, but it wasn’t.

My questions are,

  1. What is the unit of velocity? (meters? degrees? or something else?)
  2. How can I get or convert the location of the points to geographic coordinates?

My domain is 3D chunk, 105-150E(longitude), 15-50N(latitude), and 0-660km (depth).

First one is an example of my locations of points in the VTK result file.
Second one is a part of the saved test result using paraview-save data menu.
(I tested it with very low viscosity material)

screen_shot
screen_shot1920Ă—1518 215 KB

2

Kyeongjun:
ASPECT always uses the SI system of units, so coordinates are provided in meters and velocities in m/s – unless you set the “Use years in output instead of seconds” parameter to true in the input file, in which case you get m/year for velocities.

All of this also applies if you work with the “chunk” or “spherical shell” or similar geometries. There are specific in terms of spherical coordinates, but they are still only considered part of the x-y-z space. ASPECT internally doesn’t treat these domains any differently than any other domain, and so coordinates are still x-y-z. If you need the coordinates of a point in any other coordinate system (e.g. r-phi-theta), you have to do the conversion yourself, but that is of course easy with the usual set of formulas.

Best
W.

1 Like
3

Thank you for your answer.

I have one more question.

Could you explain what “Points:0, Points:1, Points:2” means?

I think it means x-y-z values in cartesian coordinates.
I can convert it to polar coordinates (r-phi-theta), but where is the reference point?

If I know where is reference point in geographical coordinates (e.g., 0,0,0 means 0°E,0°N,0m), I can convert from x-y-z values to longitude-latitude-depth values.

As you know, I want to try to interpret its result in the real Earth coordinate system.

Thank you again in advance.

Best

4

Hi Kyeongjun,

Could you explain what “Points:0, Points:1, Points:2” means?
I think it means x-y-z values in cartesian coordinates.

Yes, that refers to the x, y and z values in the cartesian system

I can convert it to polar coordinates (r-phi-theta), but where is the reference point?

If I know where is reference point in geographical coordinates (e.g., 0,0,0 means 0°E,0°N,0m), I can convert from x-y-z values to longitude-latitude-depth values.

As you know, I want to try to interpret its result in the real Earth coordinate system.

Thank you again in advance.

For reference, here is a utility function in the ASPECT that converts cartesian to spherical coordinates. That file also contains other similar functions for ellipsoidal chunks and going from spherical to cartesian coordinates.

The issue you are running into with reference frames in the spherical models (full sphere or chunk) is one we recently discussed as group. For the spherical chunk, one solution will be to allow the user to provide a reference frame (i.e., specific sets of lat, lon values) to define the chunk geometry, rather than simply lat, lon, and radius values that are not georeferenced. The ellipsoidal chunk actually currently does this.

For your models, the georeference frame will come from your input data. In other words, if the data you input into the simulation (initial conditions) is based on geographic data, you can use that rotate the spherical coordinates derived from the conversion outlined above to something that is in a georeference frame (i.e, your input data defines what the lat and lon of the corners should be, etc).

Cheers,
John

1 Like
5

Hi Kyeongjun:
to complement John’s answer: The reference point (x=0, y=0, z=0) is the center of the Earth. That should make it easy to convert into r-phi-theta.
Best
W.

1 Like
6

Hello, all.

I will try to convert it following your answers.

Thank you for your kind answers.

It will be really helpful to me.

Best,
Kyeongjun.