nav[aria-label="Primary Navigation"] { padding: 0; & ul { list-style: none; width: 100%; display: flex; flex-direction: row; justify-content: start; align-items: start; gap: 30px; padding: 0; & li { margin: 0; } & ul li { list-style: none; } } }

brainteaser viewshed wind turbines - the more you see the worse it gets...

Anonymous User

Original User: firefry82

working with arcgis 10.1 Sp1

Hi everybody,
I have a little brainteaser for you: I am trying to make a visibility analysis of wind turbines. But instead of the simple viewshed thing, I am trying to include stuff like partial visibility and multiple visibilities. The latter one is driving me insane and I hope some of you guys (or girls of course) maybe have an idea how to approach this.
I have zones around each single turbine, which indicate the level of visual distress in relation to the distance.
[ATTACH=CONFIG]26444[/ATTACH]
I can now merge them to a big viewshed which works fine.

To include the multiple visibilities I need to consider the folowing idea from "Paul et al. 2004"

[ATTACH=CONFIG]26442[/ATTACH]
Paul et al 2004 & ECOGIS
Even though this is for tension towers I think it is about the same for wind turbines. It says that the visibility increases when more when multiple turbines are visible. So he set up a factor which is 1 for the closest turbine, 1/2 for the second closest, 1/3 for the third closest and so on (1/i). Now how do I find out which turbine is closest to each cell from my raster? and then how do I get a hold of the zone-id?

This is the idea sketch:

[ATTACH=CONFIG]26447[/ATTACH]

What I tried so far:
1) converted viewsheds to polygones, made union analysis to try to intersect each one with the other to find overlapping areas. Then could go through each row and find the biggest zone-id which would be taken in account by the factor 1, the second biggest zone-id with 1/2 and so on. --> failed because union analysis can'T handle the data load

2) tried intersect the the polygon-viewsheds without the union, but then only a small part in which they ALL intersect is produced. I the searched for a intersect method that maybe finds and visualizes all the intersections from the overlaps...but I didn't find anything

3) tried to work it out with the raster calculator, but then I wouldn't know how to tell him that he should only take the zone-value when it is the closest to the cell...

I just can't get my head around it...

I am open for any ideas, phyton-scripts, ...whatever...

Image4.jpg

Image3.jpg

Image12.jpg

Image7.jpg

share_229319268894711797x.png

DougBrowning_0-1761666013232.png

RobinCoode1_0-1761664971226.png

HubsInfoSetting.JPG

Finde mehr Beiträge mit dem Tag

Accepted answers

Anonymous User

Original User: csny490

"If something is hard you should just give up, since it probably wasn't worth doing it in the first place."

NOT!

Question: What if two towers were visible: One at 200m distance, and the other at 201m distance (shouldn't that distance rank variable be a bit more scaler?)

Okay - How about this (and BTW, this would go way faster if you can coarsen the cell size). Also, make sure the cells are all alligned by setting a snap raster. Also, be aware you might fill up your RAM depending on how many pixels are involved:

1. For each tower make a viewshed/euclidean distance raster so that the output raster name has the tower id in it "tower_12", and the raster values are the distances away from the tower.
2. Convert that raster to a point featureclass.
3. Using a search cursor, load all the point FCs into a Python dictionary. The dictionary key will be a tuple of the X,Y coordinate pair.

So maybe some thing like this (warning UNTESTED!):

pointDict = {} pointFcList = ["tower_1.shp","tower_2.shp","tower_3.shp"] for pointFc in pointfcList:    searchRows = arcpy.da.SearchCursor(pointFc, ["SHAPE@XY","DISTANCE"])    towerId = pointFc.split("_")[-1]    for searchRow in searchRows:       xyKey, distance = searchRow       if xyKey not in pointDict:          pointDict[xyKey] = [(distance, towerId)]       else:          pointDict[xyKey].append((distance, towerId)) for keyKey in pointDict:    pointDict[xyKey].sort()

So now (I think) you have everything you need to compute the score:

How many towers are visible at a given x,y coordinate? len(pointDict[xyKey])
What is the closest tower's distance to a given x,y coordinate? pointDict[xyKey][0][0]
What is the closest tower id to a given x,y coordinate? pointDict[xyKey][0][1]
How about the 2nd closest tower distance? Or the 3rd? pointDict[xyKey][1][0], pointDict[xyKey][2][0]
How about the 2nd closest tower id? Or the 3rd? pointDict[xyKey][1][1], pointDict[xyKey][2][1]

All comments

Anonymous User

Original User: firefry82

okay here is a little update:
I managed to do the following
1) calculate all the intersects from each viewshed to viewshed 1, then 2 and so on.
this gives me a load of intersecting polygons for every viewshed.
2) with the help of the near-function I can now determine for each of these (intersection-)polygons how far away it is, which when brought to order gives me (theoretically) the nearest, second nearest , third nearest and so on object in relation to my first viewshed.
3) then (after conversion to raster) I can add them up by dividung the vakue of the (intersection-)raster by the number of its appearence

This works (theoretically) fine, but I have 155 viewsheds in my first set of data and if I do this calculation with every viewshed it will take a very, very long time. For an one second calculation about 7 hours. And the calculation will exceed one second by far.

How else can I approach this?

At first I just added the rasters which gave me the seemingly right outcome, but on a closer look it was skewed.

lets say the outer zone has the value of 1 while the inner (closest to the turbine) zone has the value of 6. Then if I just add up the zoned rasters six overlapping zone one's have the same value as the original zone 6. therefore I wanted to include the division

okay, the last approach seemed to be a dead end, BUT what if I would calculate the order of distance for each cell?
like:
cell 1= turbine112, turbine 7, turbine 39 ...and so on...
then I could sum up the values like this :
cell1= ([value from turbine112] *1) + ([value from turbine7] /2) + ([value from turbine39] /3)

but how do I do this? the euclidean distance only shows how far it is to the closest turbine and doesn't give me any idea which turbine it is...
damn it this shouldn't be to hard, why can't I figure this out??!

Anonymous User

Original User: curtvprice

but how do I do this? the euclidean distance only shows how far it is to the closest turbine and doesn't give me any idea which turbine it is...

The euclidean allocation grid is used to provide that information. It can be generated with its own tool or using as an optional argument to the Euclidean distance tool.

Take a deep breath, I think you're making progress. 🙂

Yes I know, I already tried it with euclidean allocation but it still only gives me the ID to the closest turbine. But I need the ID to the closest, the second closest, the thrid closest and so on. I don't see a way, how I could manipulate the tool to these specifications...

Anonymous User

Original User: csny490

But I need the ID to the closest, the second closest, the thrid closest and so on.

How about making a euclidean allocation grid for the 1st closest, then another one for the 2nd, 3rd, etc.. Then use the combine tool to basically union all the euclidean allocation grids?

Maybe it would just be easier to make a bunch of seperate viewshed rasters and then combine them together.... then cursor through the table and work out some sort of scoring system.

first of all thank you two so such for your efforts! I'm really strugelling with this and desperatly need new ideas / input.

@ csny490

How about making a euclidean allocation grid for the 1st closest, then another one for the 2nd, 3rd, etc.. Then use the combine tool to basically union all the euclidean allocation grids?

AFAIK there is no way to calculate the euclidean allocation for anything else than the closest:

[ATTACH=CONFIG]26544[/ATTACH]
this is what my euclidean allocation grid looks like (zoomed of course). It is showing precisely the which cell belongs to which turbine but it is using the SHORTEST way to do that. So all I can derive from that is that if a person would stand in that cell he would see this turbine as the closest. But I need to determine which one comes after that, and after that and so on. Because it would bother more if there are two wind turbines. BUT: this relation is NOT linear.
a possible view could look something like this:
[ATTACH=CONFIG]26546[/ATTACH]
If it wasn't for the distance I could just add up the zones, but it makes a huge difference if you see one real close (and therefore big appearing) turbine or 6 ones that are far away and therefore small appearing.

Maybe it would just be easier to make a bunch of seperate viewshed rasters and then combine them together.... then cursor through the table and work out some sort of scoring system.

I cannot do that because if I merge them, the overlays (the parts where you see more than one turbine) are lost, and if I use the union tool it crashes. And even if it wouldn't (f.e. when I do this with polygons in small steps) the union analysis takes about 3 minutes per 10 viewsheds. so even if the calculation was over after the union (which it isn't) it would take 3 mins * 15,5 *155 = 120,125 hours
If you have questions please ask. I'm glad that you are helping.

Original User: csny490

How about this:

1. For each tower make a viewshed raster and then score/reclassify it how you want (this is represented by the 1st image in your original post). So if tower id = 23, you should have a raster that is named something like "viewshed_23". For the analysis extent (the area that encompasses all the towers, reclassify your NoData pixels to 0s. You can use something like Con(IsNull("viewshed_23"), 0, "viewshed_23") to accomplish this. You have to make use of the envr "analysis extent" and "snap rastrer" settings to make sure that the calculation is done right. So in the end, you should have a bunch of rasters that all have the same rectangular extent, do not have any NoData pixles, and have only values 0,1,2,3, etc. - the 1,2,3,... representing your viewshed proiximity zones, and the 0s representing areas that are otside of the viewshed.

2. Now that all your raster have the same rectangular extent, and all the NoData pixels are now 0s, you can use the Combine tool, and end up with a table like this:

VIEWSHED_1  VIEWSHED_2  VIEWSHED_3
0           1           1               
0           2           0
0           3           0

Then you will need to add a field called "SCORE" and run an update cursor to read the VIEWSHED_* values and evaluate them. Maybe just the field calculator will suffice for this, but probably not if you have a whole bunch of towers.

ChrisSnyder

BTW: I like the look of wind towers - I think they are cool - I just don't like the sound downwind of them!

MattSayler

Maybe something from the Spatial Statistics Tools? Seems like something in there might be able to take care of the distance issue (maybe Geographically Weighted Regression?).

Anonymous User

Original User: firefry82

@csny490
Wow! It's really scary to read this, because I have done this almost the EXACT same way! Exept I used an addition where you suggested the combine-tool. I've pretty much done it in python. here are some screenshots:

the single viewsheds
[ATTACH=CONFIG]26549[/ATTACH]

the zoned-viewsheds with zeroes and big extent (btw. I don't understand what would you use the snap extent for?)
[ATTACH=CONFIG]26550[/ATTACH]

here is the python-code I used to sum (addition) the rasters together

rasters2 = arcpy.ListRasters("Raster*")
for raster in rasters:
    print "processing raster:"+raster
    #big extent
    arcpy.env.extent = "C:/studium/00_Master/clean/clipper.shp"
    #convert nodata to zero
    out1 = Con(IsNull(raster), 0, raster)

    #sum rasters together
    if i == 0:
        out2 = out1
        i += 1
    else:
        out2 = out2 + out1
        i += 1
head, tail = ntpath.split(raster)
#save final output
out2.save("C:/studium/project/script6/output.tif")

*the head and tail are for naming purposes when this script is looped.

after this I used another con to only put the values from the just generated (big) raster, which belong to the viewshed of turbine1 and so on. I have to do this because, I used the near tool to dtermine the distances from all 155 turbines to turbine 1, then 2 and so on. I sorted them asending by distance and numbered them (for the division). The only problem is: this order of appearance (which turbine comes first and which second...) is only correct for turbine1! So I do have to calculate this for all 155 turbines (155*155). In addition I need to erase the parts from the second Raster which I already calculated for Raster1. With my resolution (cell size=5 meters) a single run takes 18 minutes (on a i7 3.5GHZ with 16GB RAM and SSD). But since you would use the same method I think it won't work faster. This is why I was looking for a different solution with that euclidean distance and all.

@msayler : I don't know these tools, let me have a look at them and get back to you asap

Image5.jpg

Image6.jpg

ChrisSnyder

How about something like:

arcpy.env.extent = "C:/studium/00_Master/clean/clipper.shp"
comboGrdList = []
rastersList = arcpy.ListRasters("Raster*")
for raster in rasterList:
    out1 = Con(IsNull(raster), 0, raster)
    out1.save("twr_" + str(raster.split("_")[-1])) #assuming that raster is something like "tower_12"
    comboGrdList.append(out1)
comboGrd = arcpy.sa.Combine(comboGrdList)

I don't really understand the tower distance order weighting you are trying to implement, but the GenerateNearTable tool will give you a matrix for all 155 ^ 2 combinations and the distances between them. A dictionary that stores the distances (generated from the near table) plus an update cursor on the combine table that accesses the dictionary should be what you need I think.

Anonymous User

Original User: msayler

I'm wondering, do you really need to know the actual next towers out? I'm wondering if you could do something like:

make a raster where each cell is a count of the towers visible

create viewsheds for each tower, make 0=not visible and 1=visible (can't remember if this is what the output looks like by default)
add up the rasters per csny490's instructions, which I believe keeps the cells lined up?

create a heat map or the like (I think that's basically what you have)
normalize the values for the heat map so that they range from 0 to 1
multiply the rasters together (again, some processing may need to be done to make sure cells are coincident)

hello again, well I tried everything for the last few days but still no good results...

the combogrid

comboGrd = arcpy.sa.Combine(comboGrdList)

sounded good but it seems that combine tool doesn't work for more than 20 files, because the algorithm writes the extra info in the raster bands and 20 seems to be the maximum.

create a heat map or the like (I think that's basically what you have)
normalize the values for the heat map so that they range from 0 to 1

I tried that (not sure if I did it right), but with the multiplication small values (e.g. 0,5) keep shrinking instead of slowly growing.

I tried to resample my data to 25mx25 m then cut out the cells where the value is 0 (so no turbine visible) and then converted raster to points (see pciture1).
[ATTACH=CONFIG]26655[/ATTACH]
picture 1

I thought I could then run the "generate near table tool" on the point shapefile and the turbines shapefile (155 points or "turbines") . It took a few minutes, but actually gnereated the table. But I figure running any calculation on
this huge table (1,9 million cells * 155 locations / turbines) will take like foreever

so I am looking for a different solution on how I could accomplish this here:
[ATTACH=CONFIG]26656[/ATTACH]
picture 2: red dots are locations, grid represents (49) cells, green lines are to show the distance. Of course this is just a small scratch. Real values should be 1,9 million cells and 155 locations...

does anybody think this is even doable or do you think its time to give up on it?

Image7.jpg

Image8.jpg

Anonymous User

Original User: csny490

"If something is hard you should just give up, since it probably wasn't worth doing it in the first place."

pointDict = {} pointFcList = ["tower_1.shp","tower_2.shp","tower_3.shp"] for pointFc in pointfcList:    searchRows = arcpy.da.SearchCursor(pointFc, ["SHAPE@XY","DISTANCE"])    towerId = pointFc.split("_")[-1]    for searchRow in searchRows:       xyKey, distance = searchRow       if xyKey not in pointDict:          pointDict[xyKey] = [(distance, towerId)]       else:          pointDict[xyKey].append((distance, towerId)) for keyKey in pointDict:    pointDict[xyKey].sort()

Hey csny490, fist of all a big thank you! For the magnificent piece of code which you (if it was nothing) pulled out of your sleeve. When I try to code something its always more like a try and error thing. And then of course for the motivation. I was about to let it all go.
So here is what I did (like you advised is)
1. Set the 0-cells in my vewsheds to nodata, then resample (for testing purposes to 90x90m, 25x25m, 5x5m)
2. create euclidean distance in the same resolutions (90, 25 ,5)
3. Use con to insert the eucl distances in my modified viewshed raster
4. converted these to point shapes (snap to raster, very important!! otherwise the lenght of the dict will always be 1 and the dict will have multiple times more points
5. ran your script with these minor modifications (typos):

for pointFc in pointfcList:

must be pointFcList

for keyKey in pointDict:
pointDict[xyKey].sort()

must be pointDict[keyKey].sort()

towerId = pointFc.split("_")[-1]

I added [:-4] to cut off the ".shp"

and I have to say it works like a charm! resolution 90x90m takes only 1:30 minutes, 25x25 only 18 minutes. Sadly the 5x5 resolution takes up a sweet time and diskspace. But I guess that would have happend with any solution.

Well I think I still owe you an answer:

Question: What if two towers were visible: One at 200m distance, and the other at 201m distance (shouldn't that distance rank variable be a bit more scaler?)

No, I think. lets assume tower 1 (turbine or whatever) is in zone 6, and tower 2 in zone 5. Then the calculation would be 6 + 5/2= 8,5. And since if you have one huge wind turbine in front of you, one more, a little further away won't make too much difference.
I thought about making my zones more as you say "scale" by using the distance but I think its not worth it.

Anonymous User

Original User: csny490

Glad it worked! I replied to your message too - but yeah, go right ahead!

Nimm an der Diskussion teil

Angemeldete Mitglieder können Beiträge verfassen, Updates folgen und mehr. Neu hier? Registriere ein kostenloses Konto.

Schnellzugriffe

Community Help
Find useful guides, FAQs, and documents to help you navigate and make the most of Esri Community.