Line data Source code
1 0 : subroutine faccumulateToGrid(grid, convFuncV, nvalue, wVal,
2 : $ scaledSupport,scaledSampling,
3 : $ off, convOrigin, cfShape, loc,
4 : $ igrdpos, sinDPA, cosDPA,
5 : $ finitePointingOffset,
6 : $ doPSFOnly,
7 : $ norm,
8 0 : $ phaseGrad,
9 : $ imNX, imNY, imNP, imNC,
10 : $ cfNX, cfNY, cfNP, cfNC,
11 : $ phNX, phNY)
12 :
13 : implicit none
14 : integer imNX, imNY, imNC, imNP,
15 : $ vNX, vNY, vNC, vNP,
16 : $ cfNX, cfNY, cfNP, cfNC,
17 : $ phNX, phNY
18 :
19 : complex grid(imNX, imNY, imNP, imNC)
20 : complex convFuncV(cfNX, cfNY, cfNP, cfNC)
21 : complex nvalue
22 : double precision wVal
23 : integer scaledSupport(2)
24 : real scaledSampling(2)
25 : double precision off(2)
26 : integer convOrigin(4), cfShape(4), loc(4), igrdpos(4)
27 : double precision sinDPA, cosDPA
28 : integer finitePointingOffset, doPSFOnly
29 : complex norm
30 : complex phaseGrad(phNX, phNY)
31 :
32 : integer l_phaseGradOriginX, l_phaseGradOriginY
33 : integer iloc(4), iCFPos(4)
34 : complex wt, cfArea
35 :
36 : complex fGetCFArea
37 : integer ix,iy
38 : integer l_igrdpos(4)
39 :
40 : data iloc/1,1,1,1/, iCFPos/1,1,1,1/
41 0 : l_igrdpos(3) = igrdpos(3)+1
42 0 : l_igrdpos(4) = igrdpos(4)+1
43 : c norm=0.0
44 0 : l_phaseGradOriginX=phNX/2 + 1
45 0 : l_phaseGradOriginY=phNY/2 + 1
46 :
47 : c Currently returns 1.0
48 : cfArea = fGetCFArea(convFuncV, wVal, scaledSupport,
49 : $ scaledSampling, off, convOrigin, cfShape,
50 0 : $ cfNX, cfNY, cfNP, cfNC)
51 :
52 : c cfArea=complex(1.0, imag(cfArea))
53 : c cfArea=1.0
54 0 : do iy=-scaledSupport(2),scaledSupport(2)
55 0 : iloc(2)=nint(scaledSampling(2)*iy+off(2)-1)
56 0 : iCFPos(2)=iloc(2)+convOrigin(2)+1
57 0 : l_igrdpos(2) = loc(2)+iy+1
58 0 : do ix=-scaledSupport(1),scaledSupport(1)
59 0 : iloc(1)=nint(scaledSampling(1)*ix+off(1)-1)
60 0 : iCFPos(1) = iloc(1) + convOrigin(1) + 1
61 0 : l_igrdpos(1) = loc(1) + ix + 1
62 :
63 : wt = convFuncV(iCFPos(1), iCFPos(2),
64 0 : $ iCFPos(3), iCFPos(4))/cfArea
65 0 : if (wVal > 0.0) wt = conjg(wt)
66 :
67 0 : norm = norm + (wt)
68 :
69 : c$$$ if ((finitePointingOffset .eq. 1) .and.
70 : c$$$ $ (doPSFOnly .eq. 0)) then
71 0 : if (finitePointingOffset .eq. 1) then
72 : wt = wt * phaseGrad(iloc(1) + l_phaseGradOriginX,
73 0 : $ iloc(2) + l_phaseGradOriginY)
74 : endif
75 :
76 : grid(l_igrdpos(1), l_igrdpos(2), l_igrdpos(3), l_igrdpos(4))
77 : $ =grid(l_igrdpos(1), l_igrdpos(2),
78 0 : $ l_igrdpos(3), l_igrdpos(4)) + nvalue * wt
79 : enddo
80 : enddo
81 0 : end
82 :
83 0 : subroutine dInnerXLoop(grid,
84 : $ imNX, imNY, imNP, imNC,
85 : $ nvalue,
86 : $ scaledSupport,scaledSampling,
87 : $ convOrigin,
88 : $ iloc, l_igrdpos,iCFPos,loc,off,
89 : $ wVal,
90 : $ cfNX, cfNY,cfNP, cfNC,
91 0 : $ convFuncV,
92 : $ cfArea,
93 : $ phNX, phNY,
94 0 : $ phaseGrad,
95 : $ finitePointingOffset, doPSFOnly,
96 : $ norm)
97 : implicit none
98 : integer imNX, imNY, imNP, imNC
99 : double complex grid(imNX, imNY, imNP, imNC)
100 : complex nvalue
101 : integer scaledSupport(2), convOrigin(4)
102 : integer iloc(4), iCFPos(4),loc(4),l_igrdpos(4)
103 : double precision off(2),wVal
104 : integer cfNX, cfNY, cfNP, cfNC
105 : complex convFuncV(cfNX, cfNY, cfNP, cfNC), cfArea
106 : complex wt,norm
107 : integer finitePointingOffset, doPSFOnly
108 : integer phNX, phNY
109 : complex phaseGrad(phNX, phNY)
110 : real scaledSampling(2)
111 :
112 : integer l_phaseGradOriginX, l_phaseGradOriginY
113 :
114 : integer ix
115 :
116 0 : do ix=-scaledSupport(1),scaledSupport(1)
117 0 : iloc(1)=nint(scaledSampling(1)*ix+off(1)-1)
118 0 : iCFPos(1) = iloc(1) + convOrigin(1) + 1
119 0 : l_igrdpos(1) = loc(1) + ix + 1
120 :
121 : wt = convFuncV(iCFPos(1), iCFPos(2),
122 0 : $ iCFPos(3), iCFPos(4))/cfArea
123 0 : if (wVal > 0.0) wt = conjg(wt)
124 :
125 0 : norm = norm + (wt)
126 :
127 : c$$$ if ((finitePointingOffset .eq. 1) .and.
128 : c$$$ $ (doPSFOnly .eq. 0)) then
129 0 : if (finitePointingOffset .eq. 1) then
130 : wt = wt * phaseGrad(iloc(1) + l_phaseGradOriginX,
131 0 : $ iloc(2) + l_phaseGradOriginY)
132 : endif
133 :
134 : grid(l_igrdpos(1), l_igrdpos(2), l_igrdpos(3), l_igrdpos(4))
135 : $ =grid(l_igrdpos(1), l_igrdpos(2),
136 0 : $ l_igrdpos(3), l_igrdpos(4)) + nvalue * wt
137 : enddo
138 0 : return
139 : end
140 :
141 0 : subroutine dfaccumulateToGrid(grid, convFuncV, nvalue, wVal,
142 : $ scaledSupport,scaledSampling,
143 : $ off, convOrigin, cfShape, loc,
144 : $ igrdpos, sinDPA, cosDPA,
145 : $ finitePointingOffset,
146 : $ doPSFOnly,
147 : $ norm,
148 0 : $ phaseGrad,
149 : $ imNX, imNY, imNP, imNC,
150 : $ cfNX, cfNY, cfNP, cfNC,
151 : $ phNX, phNY)
152 :
153 : implicit none
154 : integer imNX, imNY, imNC, imNP,
155 : $ vNX, vNY, vNC, vNP,
156 : $ cfNX, cfNY, cfNP, cfNC,
157 : $ phNX, phNY
158 :
159 : double complex grid(imNX, imNY, imNP, imNC)
160 : complex convFuncV(cfNX, cfNY, cfNP, cfNC)
161 : complex nvalue
162 : double precision wVal
163 : integer scaledSupport(2)
164 : real scaledSampling(2)
165 : double precision off(2)
166 : integer convOrigin(4), cfShape(4), loc(4), igrdpos(4)
167 : double precision sinDPA, cosDPA
168 : integer finitePointingOffset, doPSFOnly
169 : complex norm
170 : complex phaseGrad(phNX, phNY)
171 :
172 : integer l_phaseGradOriginX, l_phaseGradOriginY
173 : integer iloc(4), iCFPos(4)
174 : complex wt, cfArea
175 :
176 : complex fGetCFArea
177 : integer ix,iy
178 : integer l_igrdpos(4)
179 :
180 : data iloc/1,1,1,1/, iCFPos/1,1,1,1/
181 0 : l_igrdpos(3) = igrdpos(3)+1
182 0 : l_igrdpos(4) = igrdpos(4)+1
183 : c norm=0.0
184 0 : l_phaseGradOriginX=phNX/2 + 1
185 0 : l_phaseGradOriginY=phNY/2 + 1
186 :
187 : c Currently return 1.0
188 : cfArea = fGetCFArea(convFuncV, wVal, scaledSupport,
189 : $ scaledSampling, off, convOrigin, cfShape,
190 0 : $ cfNX, cfNY, cfNP, cfNC)
191 :
192 : c cfArea=1.0
193 : c$$$!$OMP PARALLEL
194 : c$$$!$OMP& DEFAULT(NONE)
195 : c$$$!$OMP& FIRSTPRIVATE(grid)
196 : c$$$!$OMP& FIRSTPRIVATE(iloc,iCFPos,l_igrdpos,convFuncV)
197 : c$$$!$OMP& FIRSTPRIVATE(imNX,imNY,imNP,imNC,nvalue)
198 : c$$$!$OMP& FIRSTPRIVATE(scaledSupport, scaledSampling)
199 : c$$$!$OMP& FIRSTPRIVATE(convOrigin)
200 : c$$$!$OMP& FIRSTPRIVATE(loc,off,wVal,cfNX)
201 : c$$$!$OMP& FIRSTPRIVATE(cfNY, cfNP, cfNC)
202 : c$$$!$OMP& FIRSTPRIVATE(cfArea,phNX, phNY, phaseGrad)
203 : c$$$!$OMP& FIRSTPRIVATE(finitePointingOffset, doPSFOnly)
204 : c$$$!$OMP& FIRSTPRIVATE(norm)
205 : c$$$!$OMP& PRIVATE(iy) NUM_THREADS(3)
206 : c$$$!$OMP DO
207 0 : do iy=-scaledSupport(2),scaledSupport(2)
208 0 : iloc(2)=nint(scaledSampling(2)*iy+off(2)-1)
209 0 : iCFPos(2)=iloc(2)+convOrigin(2)+1
210 0 : l_igrdpos(2) = loc(2)+iy+1
211 : c$$$ call dInnerXLoop(grid, imNX, imNY, imNP, imNC,
212 : c$$$ $ nvalue, scaledSupport, scaledSampling,
213 : c$$$ $ convOrigin, iloc, l_igrdpos, iCFPos, loc,off,
214 : c$$$ $ wVal, cfNX, cfNY, cfNP, cfNC, convFuncV,cfArea,
215 : c$$$ $ phNX, phNY, phaseGrad,
216 : c$$$ $ finitePointingOffset, doPSFOnly,norm)
217 0 : do ix=-scaledSupport(1),scaledSupport(1)
218 0 : iloc(1)=nint(scaledSampling(1)*ix+off(1)-1)
219 0 : iCFPos(1) = iloc(1) + convOrigin(1) + 1
220 0 : l_igrdpos(1) = loc(1) + ix + 1
221 :
222 : wt = convFuncV(iCFPos(1), iCFPos(2),
223 0 : $ iCFPos(3), iCFPos(4))/cfArea
224 0 : if (wVal > 0.0) wt = conjg(wt)
225 :
226 : c$$$ if ((iCFPos(1) .gt. cfShape(1)) .or.
227 : c$$$ $ (iCFPos(2) .gt. cfShape(2)) .or.
228 : c$$$ $ (iCFPos(3) .gt. cfShape(3)) .or.
229 : c$$$ $ (iCFPos(4) .gt. cfShape(4))) then
230 : c$$$ write(*,*) 'O: ',iCFPos, cfShape,scaledSupport,
231 : c$$$ $ off,scaledSampling
232 : c$$$ endif
233 : c$$$ if ((iCFPos(1) .le. 0) .or.
234 : c$$$ $ (iCFPos(2) .le. 0) .or.
235 : c$$$ $ (iCFPos(3) .le. 0) .or.
236 : c$$$ $ (iCFPos(4) .le. 0)) then
237 : c$$$ write(*,*) 'U: ',iCFPos, cfShape,scaledSupport,
238 : c$$$ $ off,scaledSampling
239 : c$$$ endif
240 :
241 0 : norm = norm + (wt)
242 :
243 : c$$$ if ((finitePointingOffset .eq. 1) .and.
244 : c$$$ $ (doPSFOnly .eq. 0)) then
245 0 : if (finitePointingOffset .eq. 1) then
246 : wt = wt * phaseGrad(iloc(1) + l_phaseGradOriginX,
247 0 : $ iloc(2) + l_phaseGradOriginY)
248 : endif
249 :
250 : c$$$ if (isnan(abs(wt))) then
251 : c$$$ write(*,*) iCFPos(1), iCFPos(2), iCFPos(3), iCFPos(4),
252 : c$$$ $ cfArea,norm
253 : c$$$ return
254 : c$$$ endif
255 :
256 :
257 :
258 : grid(l_igrdpos(1), l_igrdpos(2), l_igrdpos(3), l_igrdpos(4))
259 : $ =grid(l_igrdpos(1), l_igrdpos(2),
260 0 : $ l_igrdpos(3), l_igrdpos(4)) + nvalue * wt
261 :
262 : enddo
263 : enddo
264 : c$$$!$OMP END DO
265 : c$$$!$OMP END PARALLEL
266 0 : end
267 :
268 0 : complex function fGetCFArea(convFuncV, wVal, support, sampling,
269 : $ off, convOrigin, cfShape,
270 : $ cfNX, cfNY, cfNP, cfNC)
271 :
272 : implicit none
273 : integer cfNX, cfNY, cfNP, cfNC
274 : complex convFuncV(cfNX, cfNY, cfNP, cfNC)
275 : double precision wVal
276 : integer support(2)
277 : real sampling(2)
278 : double precision off(2)
279 : integer convOrigin(4), cfShape(4)
280 :
281 : complex area, wt
282 : integer iloc(4), iCFPos(4)
283 : integer ix,iy
284 : data iloc/1,1,1,1/, iCFPos/1,1,1,1/
285 :
286 0 : fGetCFArea = 1.0
287 0 : return
288 : area=0.0
289 :
290 : do iy=-support(2), support(2)
291 : iloc(2) = nint((sampling(2)*iy+off(2))-1)
292 : iCFPos(2) = iloc(2) + convOrigin(2) + 1
293 : do ix=-support(1), support(2)
294 : iloc(1) = nint((sampling(1)*ix+off(1))-1)
295 : iCFPos(1) = iloc(1) + convOrigin(1) + 1
296 :
297 : wt = convFuncV(iCFPos(1), iCFPos(2), iCFPos(3), iCFPos(4))
298 : if (wVal > 0.0) wt = conjg(wt)
299 : area = area + wt
300 : enddo
301 : enddo
302 :
303 : fGetCFArea = area
304 : return
305 : end
306 :
|
