Problem with z-axis orientation iterating through regions on a ~200km-long chunk

[andyroo]

I have a problem that I think is related to the chunk transform matrix  orientation being tangential to a part of a 200km long fairly linear survey (elliptical hghts). I am iterating through a shapefile setting the region to shapevile bounds with an arbitrary z range (thanks to other posts on region manipulation in the forum).

I can successfully set XY to GCS but Z-axis seems like it's tangential to the overall pointcloud, so my z bounds work for some parts of the survey (generally in the middle) but end up above the ground plane at the end, and I can't figure out how to refine z orientation and properly (un)project on a per-region basis. Here is my code (although mostly it's code from other people I hacked together, and I think much of it contributed by Alexey in some form, thanks again).

Code: [Select]

import Metashape
import math
from pathlib import Path

#Set script variables below

#INPUT SHAPEFILE
in_shp = "input.shp"

#interval to step through shapefile (every nth shape) and produce dense clouds
shpStep = 20

#minimum elevation of region (meters)
zmin = -50

#max elevation of region (meters)
zmax = 1000

#initialize Metashape variables
doc = Metashape.app.document
chunk = doc.chunk
region = chunk.region
T = chunk.transform.matrix
S = chunk.transform.scale
m = Metashape.Vector([10E+10, 10E+10, 10E+10])
M = -m

#First rotate the region to the dataframe coordinate system
v_t = T * Metashape.Vector( [0,0,0,1] )
v_t.size = 3
R = chunk.crs.localframe(v_t) * T
region.rot = R.rotation().t()
chunk.region = region

#import shapes
chunk.importShapes(in_shp)

#Now iterate through selected shapes, reshape region, make dense cloud, export dense cloud
#before we get started, define a z range list from zmin and zmax
zrange = [zmin, zmax]

for shape in chunk.shapes.shapes[::shpStep]:
        print('changing region to',shape.attributes['Name'])
        #(re)initialize m and -m to really big and really small values for size and center calcs for each poly
        m = Metashape.Vector([10E+10, 10E+10, 10E+10])
        M = -m
        for vertex in shape.vertices:
            #add z vals to xy from 2D shapefile and iterate through those too
            for z in zrange:
                coord = vertex
                coord.size = 3
                coord[2] = z
                coord = chunk.crs.unproject(coord)
                print(coord)
                for i in range(3):
                    m[i] = min(m[i], coord[i])
                    M[i] = max(M[i], coord[i])

        #calculate center and size
        center = (M + m) / 2
        size = M - m

        #Apply to the region
        region.center = T.inv().mulp(center)
        region.size = size * (1/S)
        #chunk.region = region
        print('region changed to',shape.attributes['Name'])
I think my problem is related to the region rotation being relative to the chunk, so the z axis is tangential to some part of the survey, but when I try the following it doesn't seem to work either:

Code: [Select]

R = chunk.crs.localframe(center) * T
region.rot = R.rotation().t()
chunk.region = region
Any advice is appreciated.

[andyroo]

still hoping for some insight on this...