Геоинформационные системы (ГИС) и Дистанционное зондирование Земли
https://gis-lab.info/forum/
...ещё одна хорошая ссылка по basemap: http://basemaptutorial.readthedocs.org/ ... index.htmlAriki писал(а):В документации всё описано. Для этих проекций углы прямоугольника явно не задаются, но можно задать width и height в единицах проекции (метрах). Если вам известны координаты углов в географических координатах, то вы можете воспользоваться pyproj, чтобы получить их в системе координат проекции.
да, всё сработало, Ariki!Ariki писал(а):Может, очистить рисунок с помощью clf() и повторить с отрисовки basemap?
Код: Выделить всё
import os
from mpl_toolkits.basemap import Basemap
import numpy as np
import matplotlib.pyplot as plt
#Disk, directories
work_disk = 'E:'
work_dir = work_disk + os.sep + '2015_For' + os.sep + '2015_for_N' + os.sep
#Change Directory
os.chdir(work_dir)
#File with points
file_ivp = "all_points_with_dates.txt"
full_path_file_ivp = work_dir + file_ivp
#Open File with points
f = open(full_path_file_ivp, 'r')
#Loop by Lines of File with points
for year_loop in range(1900,2000):
yes_or_no_points_for_this_year=False
#a pointer to the beginning of the file
f.seek(0)
#loop by words at line
nr_line=0
for line in f:
#list
list_line = line.split()
#station id
nbd_list = list_line[0]
#latitude
lat_list = list_line[1]
#longitude
lon_list = list_line[2]
#year
year_list = list_line[3]
#month
month_list = list_line[4]
if int(year_list)==year_loop:
if yes_or_no_points_for_this_year==False:
yes_or_no_points_for_this_year=True
#************************************************************************************
#Basemap, coast, continets, rivers
mp = Basemap(projection='nplaea', boundinglat=68, lon_0=130, resolution='l', area_thresh=1000)
mp.drawcoastlines(linewidth=0.2)
mp.fillcontinents(color='ivory',lake_color='aqua')
mp.drawrivers(color='aqua')
#Meridians
meridians = np.arange(0.,351.,10.)
mp.drawmeridians(meridians,labels=[1,1,1,1],linewidth=0.125,latmax=89,fontsize=8)
meridians_thin = np.arange(0.,351.,5.)
mp.drawmeridians(meridians_thin,labels=[0,0,0,0],linewidth=0.07,latmax=89)
meridians_bold = np.arange(0.,351.,90.)
mp.drawmeridians(meridians_bold,labels=[0,0,0,0],linewidth=0.25,latmax=89,dashes=[0.01,0.01])
#Parallels
parallels = np.arange(0.,90.,1.)
mp.drawparallels(parallels,linewidth=0.125)
parallels_thin = np.arange(0.,90.,0.5)
mp.drawparallels(parallels_thin,linewidth=0.07)
parallels_bold = np.arange(0.,90.,5.)
mp.drawparallels(parallels_bold,linewidth=0.25,dashes=[0.01,0.01])
#ETOPO
mp.etopo(alpha=0.37)
ax = plt.gca()
#Parallels Labels by annotation
for lat_an in range(69,89):
x_crd=310.
y_crd=float(lat_an)
txt_an=str(lat_an)+'$^\circ$'+'N'
ax.annotate(txt_an, xy=mp(x_crd, y_crd), xycoords='data', clip_on=False, size=8, ha='center', va='center')
#************************************************************************************
first_line=nr_line
print "Year:", year_loop
nr_for_this_year=0
#Year annotation
txt_year_an=str(year_loop)
#PNG file Name
png_file=work_dir+"Png_"+txt_year_an+".png"
print png_file
x, y = mp(90,63)
plt.text(x, y, txt_year_an, fontsize=12, fontweight='bold',
ha='left',va='center',color='k',
bbox=dict(facecolor='b', alpha=0.2))
#calculate plot coordinates for point
x,y = mp(float(lon_list), float(lat_list))
mp.plot(x, y, marker=".", color="blue", markersize=2, alpha=0.5)
nr_for_this_year=nr_for_this_year+1
nr_line=nr_line+1
plt.savefig(png_file)
plt.clf()
f.close()
Да, конечно, понятно... Python-пример понятен. Премного благодарен!Ariki писал(а):Экспортировать сразу в растр с геопривязкой вряд ли получится, matplotlib для этого не предназначен. Но, зная координаты углов прямоугольника карты, можно выполнить геопривязку с помощью GDAL. Вот пример, как это сделать из Python. Ещё можно использовать утилиту gdal_translate с параметрами a_ullr и a_srs.
Ariki - уважаемый, а не дадите ли ссылочку как ещё и contour map по имеющемуся grid-файлу водрузить на созданную в этой теме подложку? Если таковой опыт у Вас имеется?Ariki писал(а):Экспортировать сразу в растр с геопривязкой вряд ли получится, matplotlib для этого не предназначен. Но, зная координаты углов прямоугольника карты, можно выполнить геопривязку с помощью GDAL. Вот пример, как это сделать из Python. Ещё можно использовать утилиту gdal_translate с параметрами a_ullr и a_srs.
это я что-то не то - всё в Surfer-терминах тянет "жаргонить" - написал, извините... надо построить изолинии по регулярному растру... но я что-то уже нашёл...Ariki писал(а):К сожалению, такого опыта нет - к стыду своему, даже не знаю, что такое contour map и grid-файл.
да, конечно... спасибо!Ariki писал(а):Ну, собственно, в документации у них примеров хватает.
Ariki, приветствую Вас. А Вы в matplolib'e эллипсы на карте не рисовали? Не tissot'ы - эллипсы искажений, а надо водрузить эллипсы результатов расчётов статистических характеристик К примеру, есть такие параметры: Широта: 69, Долгота -3, магнитуда 12, угол 240, ориентация эллипса 16, оси: 24/10, угол между главной осью эллипса и вектором мат ожидания 30 градусов. Если рисовали, то какими средствами? Надо получить подобное следующему:www.pmel.noaa.gov/maillists/tmap/ferret ... 00315.htmlAriki писал(а):Может, очистить рисунок с помощью clf() и повторить с отрисовки basemap?
Ariki, спасибо Вам! Всё годится!Ariki писал(а):nickleb, вы если хотите обратиться лично, пожалуйста, кидайте ссылку на топик в личку, а то я могу прозевать сообщение. К сожалению, mention system для призывания пользователей в топик на этом форуме не предусмотрена - я, во всяком случае, не нашёл.
Эллипсы рисовать не приходилось, я вообще далёк от статистики. Правильно я понимаю, что проблема только в развороте эллипса на заданный угол? Тогда можно сделать как здесь, только предварительно преобразовать координаты центра из географических в систему координат чертежа. Стрелочку вектора придётся изобразить отдельно с помощью pyplot.arrow().
Можете ещё здесь посмотреть.