2014-11-21 6 views
2

나는 몬트 카를로 복사 전송 코드를 사용하고 있습니다.이 전송 코드는 매체를 통해 광자를 발사 한 것을 시뮬레이트하고 랜덤 워크를 통계적으로 모델링합니다. 한 번에 한 광자를 천천히 발사하므로 벡터화하고 1000 광자를 한 번에 돌리고 싶습니다.numpy 배열의 다른 행을 다르게 나누기

나는 광자가 광학 깊이 0depth 사이의 nlayers 슬라이스로 전달되는 내 슬래브를 나누었습니다. 효과적으로, 즉 nlayers + 2 개의 영역 (nlayers에 슬래브 위의 영역과 슬래브 아래 영역)이 있음을 의미합니다. 각 단계에서 각 광자가 통과하는 레이어를 추적해야합니다.

두 개의 광자가 레이어 0에서 시작한다는 것을 이미 알고 있다고 가정 해 봅시다. 한 단계 씩 이동하여 레이어 2에서 끝나고 나머지 레이어는 레이어 6에서 끝납니다.이 배열은 배열 pastpresent으로 표시됩니다. 다음과 같다 : I는 포톤 ij 통과 여부 묘사 (nlayers + 2) 열 2 행으로 배열 traveled_through를 생성 할

[[ 0 2] 
[ 0 6]] 

(끝점을 포함). 그것은 (nlayers = 10와) 같은 것을 보일 것이다 :

[[ 1 1 1 0 0 0 0 0 0 0 0 0] 
[ 1 1 1 1 1 1 1 0 0 0 0 0]] 

내가 광자 반복 개별적으로 traveled_through의 각 행을 생성하여이 작업을 수행 할 수 있습니다,하지만 그것은 오히려 천천히, 그리고 종류의 많은 광자를 실행하는 지점을 패배 한 번에, 그렇게하지 않을 것입니다.

traveled_through = np.zeros((2, nlayers)).astype(int) 
traveled_through[ : , np.min(pastpresent, axis = 1) : np.max(pastpresent, axis = 1) + ] = 1 

생각 주어진 광자의 행에, 출발 층에서 인덱스 통해 엔딩 층을 포함한 모든로 1로 설정 될 것을 다음과 같이

는 I 어레이를 정의하려 다른 사람은 나머지 0. 그러나, 나는 다음과 같은 오류가 발생합니다 :

traveled_through[ : , np.min(pastpresent, axis = 1) : np.max(pastpresent, axis = 1) + 1 ] = 1 
IndexError: invalid slice 

내 추측은 배열의 다른 행이이 방법을 사용하여 다른 색인을 허용하지 않는 NumPy와입니다. 아무도 임의의 수의 광자와 임의의 수의 레이어에 대해 traveled_through을 생성하는 방법에 대한 제안이 있습니까?

+0

이 http://stackoverflow.com/questions/26130447/instantiate-a-matrix-with-x-zeros의 약간의 변화처럼 보이는 : 여러 행으로

, 당신은으로 작동 할 수있다 -and-the-rest-ones –

답변

2

두 광자가 항상 0에서 시작하면 다음과 같이 배열을 구성 할 수 있습니다.

먼저 설정 변수 ...

>>> pastpresent = np.array([[0, 2], [0, 6]]) 
>>> nlayers = 10 

...다음 배열을 구성 :

>>> (pastpresent[:,1][:,np.newaxis] + 1 > np.arange(nlayers+2)).astype(int) 
array([[1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0], 
     [1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0]]) 

또는 광자가 임의의 시작 층이있는 경우 :

>>> pastpresent2 = np.array([[1, 7], [3, 9]]) 
>>> (pastpresent2[:,0][:,np.newaxis] < np.arange(nlayers+2)) & 
    (pastpresent2[:,1][:,np.newaxis] + 1 > np.arange(nlayers+2)).astype(int)  
array([[0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], 
     [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0]]) 
+0

이것은 꽤 잘 작동하는 것 같습니다! 나는 마지막 코드 블록을'(pastpresent2 [:, 0] [:, np.newaxis] <= np.arange (nlayers + 2)) & (pastpresent2 [:, 1] [:, np.newaxis] + 1> np.arange (nlayers + 2)) .stype (int)'이므로 광자의 시작 레이어에 대해서도 1을 반환합니다. 감사! – DathosPachy

1

I 종류의 같은 이런 종류의 일에 대한 약간의 트릭은 logical_xor ufunc의 accumulate 방법을 포함한다 :이 b의 위치 사이의 항목을 설정하는 것이 1

>>> a = np.zeros(10, dtype=int) 
>>> b = [3, 7] 
>>> a[b] = 1 
>>> a 
array([0, 0, 0, 1, 0, 0, 0, 1, 0, 0]) 
>>> np.logical_xor.accumulate(a, out=a) 
array([0, 0, 0, 1, 1, 1, 1, 0, 0, 0]) 

에는, 먼저 인덱스 마지막 색인 배타적이므로, 정확히 무엇을했는지에 따라 1 개의 오류를 처리해야합니다.

>>> a = np.zeros((3, 10), dtype=int) 
>>> b = np.array([[1, 7], [0, 4], [3, 8]]) 
>>> b[:, 1] += 1 # handle the off by 1 error 
>>> a[np.arange(len(b))[:, None], b] = 1 
>>> a 
array([[0, 1, 0, 0, 0, 0, 0, 0, 1, 0], 
     [1, 0, 0, 0, 0, 1, 0, 0, 0, 0], 
     [0, 0, 0, 1, 0, 0, 0, 0, 0, 1]]) 
>>> np.logical_xor.accumulate(a, axis=1, out=a) 
array([[0, 1, 1, 1, 1, 1, 1, 1, 0, 0], 
     [1, 1, 1, 1, 1, 0, 0, 0, 0, 0], 
     [0, 0, 0, 1, 1, 1, 1, 1, 1, 0]])