Создать двоичную карту полноты

Я касаюсь цели моего проекта, но у меня возникает проблема: как я могу создать карту полноты? У меня есть много данных, поле с, может быть 500,000 объектов, которые представлены точками в моем участке с различным масштабом:Создать двоичную карту полноты

Without zoom First zoom Second zoom

Я хотел бы создать маску, я имею в виду, вырезать мой сюжет в крошечных пикселях и сказать, что если у меня есть объект в этом пикселе, я получаю значение: 1 (например, черный) elif, у меня нет объекта в моем пикселе, я получаю значение: 0 (например, белый).

Я создам маску, и я могу разделить каждое поле этой маской. Проблема в том, что я не знаю, как это сделать:/

Создаю первый скрипт, чтобы получить выбор по моим данным. Этот:

#!/usr/bin/python 
# coding: utf-8 

from astropy.io import fits 
from astropy.table import Table 
import numpy as np 
import matplotlib.pyplot as plt 

       ################################### 
       # Fichier contenant le champ brut # 
       ################################### 

filename = '/home/valentin/Desktop/Field52_combined_final_roughcal.fits' 

# Ouverture du fichier à l'aide d'astropy 
field = fits.open(filename) 
print "Ouverture du fichier : " + str(filename)  

# Lecture des données fits 
tbdata = field[1].data 
print "Lecture des données du fits"    

       ############################### 
       # Application du tri sur PROB # 
       ############################### 

mask = np.bitwise_and(tbdata['PROB'] < 1.1, tbdata['PROB'] > -0.1) 
new_tbdata = tbdata[mask] 
print "Création du Masque"  

      ################################################# 
      # Détermination des valeurs extremales du champ # 
      ################################################# 

# Détermination de RA_max et RA_min 
RA_max = np.max(new_tbdata['RA']) 
RA_min = np.min(new_tbdata['RA']) 
print "RA_max vaut :  " + str(RA_max) 
print "RA_min vaut :  " + str(RA_min) 

# Détermination de DEC_max et DEC_min 
DEC_max = np.max(new_tbdata['DEC']) 
DEC_min = np.min(new_tbdata['DEC']) 
print "DEC_max vaut : " + str(DEC_max) 
print "DEC_min vaut : " + str(DEC_min) 

      ######################################### 
      # Calcul de la valeur centrale du champ # 
      ######################################### 

# Détermination de RA_moyen et DEC_moyen 
RA_central = (RA_max + RA_min)/2. 
DEC_central = (DEC_max + DEC_min)/2. 

print "RA_central vaut : " + str(RA_central) 
print "DEC_central vaut : " + str(DEC_central) 

print " " 
print " ------------------------------- " 
print " " 

     ############################## 
     # Détermination de X et de Y # 
     ############################## 


# Creation du tableau 
new_col_data_X = array = (new_tbdata['RA'] - RA_central) * np.cos(DEC_central) 
new_col_data_Y = array = new_tbdata['DEC'] - DEC_central 
print 'Création du tableau' 


# Creation des nouvelles colonnes 
col_X = fits.Column(name='X', format='D', array=new_col_data_X) 
col_Y = fits.Column(name='Y', format='D', array=new_col_data_Y) 
print 'Création des nouvelles colonnes X et Y' 


# Creation de la nouvelle table 
tbdata_final = fits.BinTableHDU.from_columns(new_tbdata.columns + col_X + col_Y) 

# Ecriture du fichier de sortie .fits 
tbdata_final.writeto('{}_{}'.format(filename,'mask')) 
print 'Ecriture du nouveau fichier mask' 

field.close()

Хорошо, это работает! Но сейчас вторая часть - это тот момент:

################################################### 
################################################### 
################################################### 

filename = '/home/valentin/Desktop/Field52_combined_final_roughcal.fits_mask' 

print 'Fichier en cours de traitement' + str(filename) + '\n' 

# Ouverture du fichier à l'aide d'astropy 
field = fits.open(filename)   

# Lecture des données fits 
tbdata = field[1].data  

figure = plt.figure(1) 
plt.plot (tbdata['X'], tbdata['Y'], '.') 
plt.show()

Вы хоть представляете, как процесс? Как я могу вырезать свой сюжет в крошечном ящике?

Спасибо! UPDATE:

После ответа от armatita, я обновил свой скрипт:

################################################### 
################################################### 
################################################### 

filename = '/home/valentin/Desktop/Field52_combined_final_roughcal.fits_mask' 

print 'Fichier en cours de traitement' + str(filename) + '\n' 

# Opening file with astropy 
field = fits.open(filename)   

# fits data reading 
tbdata = field[1].data  

##### BUILDING A GRID FOR THE DATA ######## 
nodesx,nodesy = 360,360 # PIXELS IN X, PIXELS IN Y 
firstx,firsty = np.min(tbdata['X']),np.min(tbdata['Y']) 
sizex = (np.max(tbdata['X'])-np.min(tbdata['X']))/nodesx 
sizey = (np.max(tbdata['Y'])-np.min(tbdata['Y']))/nodesy 
grid = np.zeros((nodesx+1,nodesy+1),dtype='bool') # PLUS 1 TO ENSURE ALL DATA IS INSIDE GRID 

# CALCULATING GRID COORDINATES OF DATA 
indx = np.int_((tbdata['X']-firstx)/sizex) 
indy = np.int_((tbdata['Y']-firsty)/sizey) 
grid[indx,indy] = True # WHERE DATA EXISTS SET TRUE 

# PLOT MY FINAL IMAGE 
plt.imshow(grid.T,origin='lower',cmap='binary',interpolation='nearest') 
plt.show()

Я нахожу этот сюжет:

solution

Так что, когда я играю с размером бен, я может видеть более или менее пустое значение, указывающее объект или нет в моем пикселе :)

источник

2016-03-22 Deadpool

Обычно это процесс вставки y наши данные в сетку (пиксель мудрый или узловой). В следующем примере построена сетка (2D-массив) и вычисляется «координаты сетки» для выборочных данных. После того, как он имеет эти координаты сетки (которые в действительности являются ничем, кроме индексов массива), вы можете просто установить эти элементы в значение True. Проверьте следующий пример:

import numpy as np 
import matplotlib.pyplot as plt 

x = np.random.normal(0,1,1000) 
y = np.random.normal(0,1,1000) 

##### BUILDING A GRID FOR THE DATA ######## 
nodesx,nodesy = 100,100 # PIXELS IN X, PIXELS IN Y 
firstx,firsty = x.min(),y.min() 
sizex = (x.max()-x.min())/nodesx 
sizey = (y.max()-y.min())/nodesy 
grid = np.zeros((nodesx+1,nodesy+1),dtype='bool') # PLUS 1 TO ENSURE ALL DATA IS INSIDE GRID 

# CALCULATING GRID COORDINATES OF DATA 
indx = np.int_((x-firstx)/sizex) 
indy = np.int_((y-firsty)/sizey) 
grid[indx,indy] = True # WHERE DATA EXISTS SET TRUE 

# PLOT MY FINAL IMAGE 
plt.imshow(grid.T,origin='lower',cmap='binary',interpolation='nearest') 
plt.show()

, что приводит к: