Add z-value to shapefile made of multipatch data (roofsurfaces)

[Jacob_JeffB]

Hey guys, i am pretty new in 3D analysis in ArcMap/Scene (10.5.1.)

It might be pretty simple question ( and it hopefully is).

I am working with 3D buildings of my hometown to analyze roof solar potential.

My final goal is to produce slope and aspect data.

What kind of data i use:

-   109 multipatch (.gml) files (each contains GroundSurface, RoofSurface, Wallsurface)

>> i am interested in Roofsurfaces to analyze slope and aspect

-  1 shapefile (it contains the 109 roofsurfaces)

>> i produced it myself and used it for 2D analysis (for example: footprint to estimate and distinguish flat and sloping roof types)

What i am trying to do is:

- add z_value to the shapefile

- to create a raster elevation map

- and finally do slope and aspect analysis

My problem is:

- i dont know much about multipatch and gml. data

- i dont know how to create this z-value out of my existing data and add it the attribute table of the shapefile (if had the z-value i could use polygon to raster)

> Is it even possible? Can i use the shapefile? Or do i have to use the raw data (multipatch / gml.)

I hope you can understand what i am trying to do, and what kind of data i use.

Thanks in advance,

JacobB

[XanderBakker]

I just noticed that in the screenshot before I used the DWG files that you can download at the same site which can be read directly using ArcGIS (no conversion required).

[Jacob_JeffB]

Thank you very much for helps. It unforunately still does not work. It is probably an issue with the computers/servers i am working on (at my university that provides free ArcGis Software). I think its time get in touch with the administrator.

Here is what i tried for DWG files. The tool was successfully finished. But the output data seems to be unusable.

This what i tried with the GML data. Same ERROR as above.

1

2

3

[XanderBakker]

I just wrote a script to process all the DWG files, but there seems to be a problem with some of them. When I have a (partial) result I will post it back here.

[XanderBakker · Accepted solution]

Find below the code that I wrote to process the data (DWG files):

import arcpy

def main():
    import os

    # settings
    ws = r'C:\GeoNet\CityGML\DWG_LoD2'
    fgdb = r'C:\GeoNet\CityGML\gdb\Hannover.gdb'
    dsm_name = 'DSM'
    slope_name = 'Slope'
    aspect_name = 'Aspect'
    cellsize = 1.0
    sr = arcpy.SpatialReference(25832)  # ETRS_1989_UTM_Zone_32N

    # necessary environment settngs
    arcpy.env.outputCoordinateSystem= sr
    arcpy.env.workspace = ws

    # create list of DWG file in folder (ws)
    lst = arcpy.ListFiles("*.DWG")

    # check out SA license
    arcpy.CheckOutExtension("Spatial")

    # loop through dwg file
    lst_ras = []
    lst_err = []
    for dwg in lst:
        dwg_mp = os.path.join(ws, dwg, 'MultiPatch')
        print("Processing: {}".format(dwg))

        # correct extent
        ext = arcpy.Describe(dwg_mp).extent
        ext_corr = CorrectExtent(ext, cellsize)
        arcpy.env.extent = ext_corr

        # define output raster name
        ras_name = 'ras{}'.format(os.path.splitext(dwg)[0])
        out_ras = os.path.join(fgdb, ras_name)

        # create raster from multipatch and add to lost
        try:
            # some files create a error. I will skip those
            arcpy.MultipatchToRaster_conversion(dwg_mp, out_ras, cellsize)
            lst_ras.append(out_ras)
        except Exception as e:
            print "ERROR: ", e
            lst_err.append(dwg)


    print "DWG that produced Errors:"
    print "\n - ".join(lst_err)

    # mosaic raster into new raster
    print("mosaic into new raster...")
    input_ras = ';'.join(lst_ras)
    arcpy.MosaicToNewRaster_management(input_ras, fgdb, dsm_name, sr, "32_BIT_FLOAT", cellsize, 1, "MEAN")

    # set extent to output DSM
    dsm_ras = os.path.join(fgdb, dsm_name)
    arcpy.env.extent = dsm_ras

    # calculate slope
    print("calculate slope...")
    slope_ras = os.path.join(fgdb, slope_name)
    out_slope = arcpy.sa.Slope(dsm_ras, "PERCENT_RISE")  # or use DEGREE
    out_slope.save(slope_ras)

    # calculate aspect
    print("calculate aspect...")
    aspect_ras = os.path.join(fgdb, aspect_name)
    out_aspect = arcpy.sa.Aspect(dsm_ras)
    out_aspect.save(aspect_ras)


def CorrectExtent(ext, cellsize):
    ext.XMin = divmod(ext.XMin, cellsize)[0]
    ext.YMin = divmod(ext.YMin, cellsize)[0]
    ext.XMax = divmod(ext.XMax, cellsize)[0] + cellsize
    ext.YMax = divmod(ext.YMax, cellsize)[0] + cellsize
    return ext

if __name__ == '__main__':
    main()
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

The file geodatabase with the resulting rasters (DSM, slope, aspect), is about 1.5GB (compressed a little over 210MB). You can download it here: HannoverCityGML.zip - Google Drive 

[Jacob_JeffB]

Thank you so much!

[XanderBakker]

You are welcome