[GRASS-user] Issue with r.sun.daily in combination with multiprocessing

2017-05-15 Thread Tinsley Wolfs
Hi list,

I am trying to build a script to run r.sun.daily so it only runs for the center 
day of each month of a year. However, when the code starts r.sun.daily script, 
it creates 4 processes for the first day but only 2 processes for any following 
day. Even though I set ' nprocs=4' when calling the script.

I am trying to use a similar setup as the r.horizon parallel script found 
here.
 This script uses the multiprocessing module to have the script wait until the 
process is finished before feeding new parameters to a worker. I edited this so 
it changes the 'start_day' and 'end_day' parameters, but only creates one 
worker that launches the r.sun.daily script with 'nprocs=4', then waits until 
this run is done using '.join', before providing the worker with a different 
start and end day and starting r.sun.daily again with ' nprocs=4'.

What could be the issue that it switches to two processes after the 1st day?

Regards,
Tinsley
#!/usr/bin/env python
import os, sys, shutil, re
import datetime as dt
import multiprocessing as mp
import gdal, ogr, osr

gisbase = os.environ['GISBASE'] = r'C:\\Program Files\\GRASS GIS 7.2.1\\'
gisdbase = r'C:\data\grass'
gisaddon = os.environ['GRASS_ADDON_PATH'] = 
r'C:\Users\*\AppData\Roaming\GRASS7\addons\scripts\\'
sys.path.append(os.path.join(os.environ['GISBASE'], "etc", "python"))

# import grass scripts
import grass.script as grass
import grass.script.setup as gsetup

# days of the year to model in r.sun
days = []
day1 = dt.datetime(2015, 1, 15).timetuple().tm_yday
days.append(day1)
day2 = dt.datetime(2015, 2, 14).timetuple().tm_yday
days.append(day2)
day3 = dt.datetime(2015, 3, 15).timetuple().tm_yday
days.append(day3)
day4 = dt.datetime(2015, 4, 15).timetuple().tm_yday
days.append(day4)
day5 = dt.datetime(2015, 5, 15).timetuple().tm_yday
days.append(day5)
day6 = dt.datetime(2015, 6, 15).timetuple().tm_yday
days.append(day6)
day7 = dt.datetime(2015, 7, 15).timetuple().tm_yday
days.append(day7)
day8 = dt.datetime(2015, 8, 15).timetuple().tm_yday
days.append(day8)
day9 = dt.datetime(2015, 9, 15).timetuple().tm_yday
days.append(day9)
day10 = dt.datetime(2015, 10, 15).timetuple().tm_yday
days.append(day10)
day11 = dt.datetime(2015, 11, 15).timetuple().tm_yday
days.append(day11)
day12 = dt.datetime(2015, 12, 15).timetuple().tm_yday
days.append(day12)

# sample area
sampleList1 = ['lavbebygNet_024359']
clusterList = []
clusterList.append(sampleList1)


# open log file
global lf
tlog = dt.datetime.strftime(dt.datetime.now(),"%Y-%m-%d_h%Hm%M")
lf = open('rsun_'+tlog+'.log', 'w')

# parameters
location = 'Denmark'
mhorizon = 'Solar'
maxDistance = '500'
hstep = 15
distance = '1.0'

def worker_rsun(vectorId,dayStart,dayEnd):
# setup grass environment, including temporary mapset
mtemp = 'temp'+str(1).zfill(2)
gsetup.init(gisbase, gisdbase, location, mtemp)
grass.run_command("g.region",vector = vectorId +"@Solar", res='0.4')

# map names
prefsun = vectorId+'_glob'
demh = vectorId+'@Solar'
horizonBase = vectorId+'hor'
aspectIn = vectorId+'aspect'+'@Solar'
slopeIn = vectorId+'slope'+'@Solar'

# run r.horizon for one azimuth angle
dayStart = dayStart
dayEnd = dayEnd
grass.run_command("r.sun.daily", elevation=demh, \
  start_day=dayStart, end_day=dayEnd, aspect=aspectIn, \
  slope=slopeIn, horizon_basename=horizonBase, \
  horizon_step=15, day_step=1, 
glob_rad_basename=prefsun, \
  nprocs=4, overwrite=True)

def create_temp():
global temp
# create temporary matsets
gsetup.init(gisbase, gisdbase, location, 'PERMANENT')
temp = 'temp'+str(1).zfill(2)
temp_path = gisdbase+'/'+location+'/'+temp
if os.path.exists(temp_path) == False:
grass.run_command("g.mapset", flags="c", mapset = temp, quiet = 
1)
printout(temp+' mapset created.')
else:
printout(temp+' mapset already exists. skipping...')

def remove_temp():
# setup is in the target mapset, not temp mapsets
gsetup.init(gisbase, gisdbase, location, mhorizon)

# delete the temp mapset
mapset_temp = 'temp'+str(1).zfill(2)
grass.run_command("g.mapsets", mapset=mapset_temp, operation = 'remove')
temp_path = gisdbase+'/'+location+'/'+mapset_temp
shutil.rmtree(temp_path)

def copy_mapset(mapset_from,mapset_to,vectorId):
# list contents of temp mapset
grass.run_command("g.mapset", mapset=mapset_from, quiet=1)
raster_list = grass.list_strings('rast', pattern=vectorId+'_glob')
# switch to target mapset and copy rasters over
grass.run_command("g.mapset", mapset=mapset_to,quiet=1)
for rast in raster_list:
cmd = 

[GRASS-user] issues with r.horizon incombination with multiprocessing

2017-05-02 Thread Tinsley Wolfs
Hi list,

Based on a discussion in the list archive, I copied and edited a script that 
parallelizes the r.horizon module across multiple cores 
(link).
 I have updated this script to work with the current module parameter naming 
schemes where required. Additionally I added a for loop to the main() function, 
as I want it to do the r.horizon module for several areas sequentially (setting 
g.region, etc...).

I have several issues:

1. Each time it creates the workers for the last r.horizon angles for a list 
item, it also creates the workers for the first angles of the next item in the 
list. But it should finish the first workers, copy the raster from the 
temporary mapsets to the current mapset and delete the temporary mapsets before 
creating new workers and starting the process again.
2. It doesn't calculate the last angles of the last area part of the list that 
it loops over. Though I think this is because of the first issue which is that 
the order of execution is off.

Why doesn't it follow the correct order as per the script but move on to the 
next step before finishing the last? I run the script from the simple python 
editor in Grass 7.2.0.
#!/usr/bin/env python

import sys, os, shutil, re
import datetime as dt
import multiprocessing as mp
import subprocess
import atexit

gisbase = os.environ['GISBASE'] = r'C:\\Program Files\\GRASS GIS 7.2.0\\'
gisdbase = r'C:\data\grass'
sys.path.append(os.path.join(os.environ['GISBASE'], "etc", "python"))

import grass.script as grass
import grass.script.setup as gsetup

# code to create list with name of vectors to iterate over
className = "bykerne" + "Net"
mapName = className + '@Solar'  # counting layers (sample clusters)
result = grass.read_command('v.info',   # in vector file
flags = 'e',
map = mapName)

result = str.split(result, '\n')
layerCount = result[12] # dblink_num = layers 
layerCount = str.split(layerCount, '=')
layerCount = layerCount[1]  # number of layers
layerCount = int(layerCount)

layerN = 1  # start layer count
outputName = className + '_' + str(layerN - 1)
vectorList = grass.read_command('g.list', type = 'vector', mapset = 'Solar')
vectorList = str.split(vectorList, '\n')
vectorList = ''.join(vectorList)
vectorList = str.split(vectorList, '\r')
vectorList.pop(0)
vectorList.pop(0)
vectorList.pop(-1) # final vector list

# parameters
location = 'Denmark'
mhorizon = 'Solar'  
maxDistance = '500' 
hstep = 15  
distance = '1.0'

# start of function definitions:

def worker_horizon(cpu,demh,azimuth, outputNameS):

# Setup Grass Environment, including temporary mapset

mtemp = 'temp'+str(cpu).zfill(2)

gsetup.init(gisbase, gisdbase, location, mtemp)

grass.run_command("g.region",vector = outputNameS +"@Solar")



# map names

prefhor = demh+'hor'

demh = demh+'@Solar'



# run r.horizon for one azimuth angle

azi = str(azimuth)

grass.run_command("r.horizon", elevation=demh, \

  direction=azi, maxdistance=maxDistance, \

  output=prefhor, dist=distance, overwrite = True)



# rename horizon map to proper azimuth name

az = str(azimuth).zfill(3)

tempmap = prefhor+'_0'

horizonmap = prefhor+'_'+az

cmd = tempmap+","+horizonmap

grass.run_command("g.rename",rast=cmd,overwrite = True)



def create_temp(cores):

gsetup.init(gisbase, gisdbase, location, 'PERMANENT')

for count in range(0,cores,1):

temp = 'temp'+str(count).zfill(2)

temp_path = gisdbase+'/'+location+'/'+temp

if(os.path.exists(temp_path) == False):

grass.run_command("g.mapset",flags="c", mapset=temp, 
quiet=1)

print(temp+" mapset created.")

else:

print(temp+" mapset already exists. skipping...")



def remove_temp(cores):

# setup is in the target mapset, not temp mapsets

gsetup.init(gisbase, gisdbase, location, mhorizon)

# Delete the temp mapset

for count in range(0,cores):

mapset_temp = 'temp'+str(count).zfill(2)

grass.run_command("g.mapsets", mapset=mapset_temp, operation = 
'remove')

temp_path = gisdbase+'/'+location+'/'+mapset_temp

shutil.rmtree(temp_path)

if os.path.exists(temp_path):
shutil.rmtree(temp_path)


def copy_mapset(mapset_from,mapset_to,regfilter,overwrite):

# list contents of temp mapset

grass.run_command("g.mapset", mapset=mapset_from, quiet=1)

raster_list =