|
1 module realizer_mod |
|
2 |
|
3 ! Realizer module written for the MidPtFM2d program (now called realizer) by |
|
4 ! Forrest Hoffman on Fri May 4 16:56:34 EDT 2007 to eliminate duplicated |
|
5 ! parameters and common blocks. Parts of the code were converted to |
|
6 ! Fortran-90 style and intrinsics. Compiled and tested with gfortran |
|
7 ! version 4.1.1 under Fedora Core 6 Linux. |
|
8 |
|
9 implicit none |
|
10 save |
|
11 |
|
12 ! Parameters |
|
13 integer, parameter :: r8 = selected_real_kind(12) ! 8 byte real |
|
14 integer, parameter :: r4 = selected_real_kind( 6) ! 4 byte real |
|
15 integer, parameter :: nr = kind(1.0) ! native real |
|
16 integer, parameter :: i8 = selected_int_kind(13) ! 8 byte integer |
|
17 integer, parameter :: i4 = selected_int_kind( 6) ! 4 byte integer |
|
18 integer, parameter :: ni = kind(1) ! native integer |
|
19 integer(i4), parameter :: maxN=515 |
|
20 integer(i4), parameter :: maxsp=10 |
|
21 real(r8), parameter :: rad_per_deg=(3.14159265358979323846/180.0) |
|
22 real(r8), parameter :: deg_per_rad=(180.0/3.14159265358979323846) |
|
23 |
|
24 ! Global Variables |
|
25 integer(i4) :: numcells |
|
26 integer(i4) :: itiedsplit |
|
27 integer(i4) :: irowcols |
|
28 integer(i4) :: mtsplit |
|
29 integer(i4) :: mtresolved |
|
30 integer(i4) :: isupply |
|
31 integer(i4) :: ihabs |
|
32 integer(i4) :: maxgray |
|
33 real(r4) :: northeast |
|
34 real(r4) :: southeast |
|
35 real(r4) :: southwest |
|
36 real(r4) :: northwest |
|
37 |
|
38 ! Private Member Functions |
|
39 public dump_map |
|
40 public mpfm2d |
|
41 private StackSort ! unused? |
|
42 private sort ! unused? |
|
43 private SwapDown |
|
44 public Interchange |
|
45 public HeapSort |
|
46 private SD2 |
|
47 private I2 |
|
48 public HS2 |
|
49 private Gauss |
|
50 private cosd |
|
51 private sind |
|
52 private atand |
|
53 private atan2d |
|
54 |
|
55 contains |
|
56 ! |
|
57 !------------------------------------------------------------------------------ |
|
58 ! |
|
59 subroutine dump_map(P) |
|
60 implicit none |
|
61 integer(i4), intent(in) :: P(maxN*maxN,0:maxsp+3) |
|
62 ! local variables |
|
63 integer(i4) :: i, j ! loop indices |
|
64 integer(i4) :: temp(maxsp+3) |
|
65 character (len=40) :: fmt |
|
66 ! |
|
67 write(fmt,'("(i5, 1x, ",I0,"(i5, 1x), 1x, 3(i5, 1x))")') ihabs+1 |
|
68 do i=1,numcells |
|
69 do j = 0,ihabs+3 |
|
70 temp(j) = P(i,j) |
|
71 end do |
|
72 print fmt, i, (temp(j), j=0,ihabs+3) |
|
73 end do |
|
74 ! |
|
75 return |
|
76 end subroutine dump_map |
|
77 ! |
|
78 !------------------------------------------------------------------------------ |
|
79 ! |
|
80 real(r4) function f3(delta, x0, x1, x2, sigma) |
|
81 implicit none |
|
82 real(r4), intent(in) :: delta |
|
83 real(r4), intent(in) :: x0, x1, x2 |
|
84 real(r4), intent(in) :: sigma |
|
85 |
|
86 f3 = (x0+x1+x2) / 3.0 + delta * Gauss(0.,sigma) |
|
87 |
|
88 return |
|
89 end function f3 |
|
90 ! |
|
91 !------------------------------------------------------------------------------ |
|
92 ! |
|
93 real(r4) function f4(delta, x0, x1, x2, x3, sigma) |
|
94 implicit none |
|
95 real(r4), intent(in) :: delta |
|
96 real(r4), intent(in) :: x0, x1, x2, x3 |
|
97 real(r4), intent(in) :: sigma |
|
98 |
|
99 f4 = (x0+x1+x2+x3) / 4.0 + delta * Gauss(0.,sigma) |
|
100 |
|
101 return |
|
102 end function f4 |
|
103 ! |
|
104 !------------------------------------------------------------------------------ |
|
105 ! |
|
106 subroutine mpfm2d(X, maxlevel, sigma, H, addition, wrap, gradient, iz) |
|
107 implicit none |
|
108 real(r4), intent(inout) :: X(0:maxN,0:maxN) |
|
109 integer(i4), intent(in) :: maxlevel |
|
110 real(r4), intent(in) :: sigma |
|
111 real(r4), intent(in) :: H(maxN,maxsp) |
|
112 logical, intent(in) :: addition |
|
113 logical, intent(in) :: wrap |
|
114 logical, intent(in) :: gradient |
|
115 integer(i4), intent(in) :: iz |
|
116 ! local variables |
|
117 integer(i4) :: N |
|
118 real(r4) :: delta |
|
119 integer(i4) :: DD |
|
120 integer(i4) :: d |
|
121 integer(i4) :: is, ix, iy ! loop indices |
|
122 ! |
|
123 print *, "Generating fractal" |
|
124 N = 2**maxlevel |
|
125 delta = sigma |
|
126 X(0,0) = delta * Gauss(0.,sigma) |
|
127 X(0,N) = delta * Gauss(0.,sigma) |
|
128 X(N,0) = delta * Gauss(0.,sigma) |
|
129 X(N,N) = delta * Gauss(0.,sigma) |
|
130 ! |
|
131 if (gradient) then |
|
132 X(0,0) = southwest |
|
133 X(0,N) = southeast |
|
134 X(N,0) = northwest |
|
135 X(N,N) = northeast |
|
136 end if |
|
137 ! |
|
138 DD = N |
|
139 d = N / 2 |
|
140 ! |
|
141 do is = 1, maxlevel |
|
142 ! going from grid type I to type II |
|
143 delta = delta * 0.5**(0.5 * H(is,iz)) |
|
144 do ix = d, N-d, DD |
|
145 do iy = d, N-d, DD |
|
146 X(ix,iy) = f4(delta, X(ix+d,iy+d), X(ix+d,iy-d), & |
|
147 X(ix-d,iy+d), X(ix-d,iy-d), sigma) |
|
148 end do |
|
149 end do |
|
150 ! displace other points also if needed |
|
151 if (addition) then |
|
152 do ix = 0, N, DD |
|
153 do iy = 0, N, DD |
|
154 X(ix,iy) = X(ix,iy) + delta * Gauss(0.,sigma) |
|
155 end do |
|
156 end do |
|
157 end if |
|
158 ! going from grid type II to type I |
|
159 delta = delta * 0.5**(0.5*H(is,iz)) |
|
160 ! interpolate and offset boundary grid points |
|
161 do ix = d, N-d, DD |
|
162 X(ix,0) = f3(delta, X(ix+d,0), X(ix-d,0), X(ix,d), sigma) |
|
163 X(ix,N) = f3(delta, X(ix+d,N), X(ix-d,N), X(ix,N-d), sigma) |
|
164 X(0,ix) = f3(delta, X(0,ix+d), X(0,ix-d), X(d,ix), sigma) |
|
165 X(N,ix) = f3(delta, X(N,ix+d), X(N,ix-d), X(N-d,ix), sigma) |
|
166 if (wrap) then |
|
167 X(ix,N) = X(ix,0) |
|
168 X(N,ix) = X(0,ix) |
|
169 end if |
|
170 end do |
|
171 ! interpolate and offset inter grid points |
|
172 do ix = d, N-d, DD |
|
173 do iy = DD, N-d, DD |
|
174 X(ix,iy) = f4(delta, X(ix,iy+d), X(ix,iy-d), & |
|
175 X(ix+d,iy), X(ix-d,iy), sigma) |
|
176 end do |
|
177 end do |
|
178 do ix = DD, N-d, DD |
|
179 do iy = d, N-d, DD |
|
180 X(ix,iy) = f4(delta, X(ix,iy+d), X(ix,iy-d), & |
|
181 X(ix+d,iy), X(ix-d,iy), sigma) |
|
182 end do |
|
183 end do |
|
184 ! displace other points also if needed |
|
185 if (addition) then |
|
186 do ix = 0, N, DD |
|
187 do iy = 0, N, DD |
|
188 X(ix,iy) = X(ix,iy) + delta * Gauss(0.,sigma) |
|
189 end do |
|
190 end do |
|
191 do ix = d, N-d, DD |
|
192 do iy = d, N-d, DD |
|
193 X(ix,iy) = X(ix,iy) + delta * Gauss(0.,sigma) |
|
194 end do |
|
195 end do |
|
196 end if |
|
197 ! |
|
198 DD = DD / 2 |
|
199 d = d / 2 |
|
200 ! |
|
201 end do ! for all levels |
|
202 return |
|
203 end subroutine mpfm2d |
|
204 ! |
|
205 !------------------------------------------------------------------------------ |
|
206 ! |
|
207 real(r4) function Gauss(xm, xs) |
|
208 implicit none |
|
209 real(r4), intent(in) :: xm |
|
210 real(r4), intent(in) :: xs |
|
211 ! local variables |
|
212 real(r4) :: total |
|
213 real(r4) :: p |
|
214 integer(i4) :: i |
|
215 ! |
|
216 ! create a normal number by summing 12 uniform random numbers |
|
217 total = 0.0 |
|
218 do i = 1,12 |
|
219 call random_number(p) |
|
220 total = total + p |
|
221 end do |
|
222 Gauss = (total - 6.0) * xs + xm |
|
223 ! |
|
224 return |
|
225 end function Gauss |
|
226 ! |
|
227 !------------------------------------------------------------------------------ |
|
228 ! |
|
229 subroutine StackSort(n,arr) |
|
230 implicit none |
|
231 integer(i4), intent(in) :: n |
|
232 real(r4), intent(inout) :: arr(n) |
|
233 ! local variables |
|
234 integer(i4) :: i, j ! loop indices |
|
235 real(r4) :: a |
|
236 ! |
|
237 ! pick out each element in turn |
|
238 do j = 2,n |
|
239 a = arr(j) |
|
240 ! look for the place to insert it |
|
241 do i = j-1,1,-1 |
|
242 if (arr(i) <= a) go to 10 |
|
243 arr(i+1) = arr(i) |
|
244 end do |
|
245 i = 0 |
|
246 ! insert the number here |
|
247 10 arr(i+1) = a |
|
248 end do |
|
249 return |
|
250 end subroutine StackSort |
|
251 ! |
|
252 !------------------------------------------------------------------------------ |
|
253 ! |
|
254 ! sort routine from Numerical Recipes |
|
255 subroutine sort (n, ra) |
|
256 implicit none |
|
257 integer(i4), intent(in) :: n |
|
258 real(r4), intent(inout) :: ra(n) |
|
259 ! local variables |
|
260 integer(i4) :: i, j, l |
|
261 integer(i4) :: ir |
|
262 real(r4) :: rra |
|
263 |
|
264 l = n/2 + 1 |
|
265 ir = n |
|
266 ! |
|
267 10 continue |
|
268 if (l .gt. 1) then |
|
269 l = l-1 |
|
270 rra = ra(l) |
|
271 else |
|
272 rra = ra(ir) |
|
273 ra(ir) = ra(1) |
|
274 ir = ir-1 |
|
275 if (ir .eq. 1) then |
|
276 ra(1) = rra |
|
277 return |
|
278 end if |
|
279 end if |
|
280 i = l |
|
281 j = l+l |
|
282 ! |
|
283 20 if (j .le. ir) then |
|
284 if (j .lt. ir) then |
|
285 if (ra(j) .lt. ra(j+1)) j = j+1 |
|
286 end if |
|
287 if (rra .lt. ra(j)) then |
|
288 ra(i) = ra(j) |
|
289 i = j |
|
290 j = j+j |
|
291 else |
|
292 j = ir+1 |
|
293 end if |
|
294 goto 20 |
|
295 end if |
|
296 ! |
|
297 ra(i) = rra |
|
298 goto 10 |
|
299 ! |
|
300 end subroutine sort |
|
301 ! |
|
302 !------------------------------------------------------------------------------ |
|
303 ! |
|
304 subroutine SwapDown(Heap, r, isize, istart, keyfield) |
|
305 implicit none |
|
306 integer(i4), intent(inout) :: Heap(maxN*maxN,0:maxsp+3) |
|
307 integer(i4), intent(inout) :: r |
|
308 integer(i4), intent(in) :: isize |
|
309 integer(i4), intent(in) :: istart |
|
310 integer(i4), intent(in) :: keyfield |
|
311 ! local variables |
|
312 integer(i4) :: child |
|
313 integer(i4) :: ipos |
|
314 logical :: done |
|
315 ! |
|
316 ! print *, 'entering SwapDown' |
|
317 ! print *, 'params are ', r, isize, istart, keyfield |
|
318 ipos = istart - 1 |
|
319 done = .false. |
|
320 child = 2 * r |
|
321 do while (.not. done .and. child .le. isize) |
|
322 ! find the largest child |
|
323 if (child .lt. isize .and. Heap(child+ipos,keyfield) .lt. & |
|
324 Heap(child+ipos+1,keyfield)) child = child + 1 |
|
325 ! interchange node and largest child if necessary |
|
326 ! and move down to the next subtree |
|
327 if (Heap(r+ipos,keyfield) .lt. Heap(child+ipos,keyfield)) then |
|
328 Call Interchange(Heap,r+ipos,child+ipos) |
|
329 r = child |
|
330 child = 2 * child |
|
331 else |
|
332 done = .true. |
|
333 end if |
|
334 ! |
|
335 end do |
|
336 ! |
|
337 return |
|
338 end subroutine SwapDown |
|
339 ! |
|
340 !------------------------------------------------------------------------------ |
|
341 ! |
|
342 subroutine Interchange(X,m,n) |
|
343 implicit none |
|
344 integer(i4), intent(inout) :: X(maxN*maxN,0:maxsp+3) |
|
345 integer(i4), intent(in) :: m |
|
346 integer(i4), intent(in) :: n |
|
347 ! local variables |
|
348 integer(i4) :: i |
|
349 integer(i4) :: temp |
|
350 ! |
|
351 do i = 0, ihabs+3 |
|
352 temp = X(m,i) |
|
353 X(m,i) = X(n,i) |
|
354 X(n,i) = temp |
|
355 end do |
|
356 ! |
|
357 return |
|
358 end subroutine Interchange |
|
359 ! |
|
360 !------------------------------------------------------------------------------ |
|
361 ! |
|
362 subroutine HeapSort(X,isize,istart,keyfield) |
|
363 implicit none |
|
364 integer(i4), intent(inout) :: X(maxN*maxN,0:maxsp+3) |
|
365 integer(i4), intent(in) :: isize |
|
366 integer(i4), intent(in) :: istart |
|
367 integer(i4), intent(in) :: keyfield |
|
368 ! local variables |
|
369 integer(i4) :: i ! loop indices |
|
370 integer(i4) :: r |
|
371 integer(i4) :: tempr |
|
372 ! |
|
373 ! print *, "isize in HeapSort is ", isize |
|
374 ! print *, "istart in HeapSort is ", istart |
|
375 ! print *, "keyfield in HeapSort is ", keyfield |
|
376 ! heapify |
|
377 do r = isize/2, 1, -1 |
|
378 tempr = r |
|
379 !Call SwapDown(X,r,isize,istart,keyfield) |
|
380 Call SwapDown(X,tempr,isize,istart,keyfield) |
|
381 !r = tempr |
|
382 end do |
|
383 ! |
|
384 ! print *, "P after heapify is " |
|
385 ! do i=istart,isize |
|
386 ! print *, i, X(i,keyfield) |
|
387 ! end do |
|
388 ! |
|
389 ! repeatedly put largest item in root at end of list, |
|
390 ! prune it from the tree |
|
391 ! and apply SwapDown to the rest of the tree |
|
392 ! |
|
393 do i = isize, 2, -1 |
|
394 Call Interchange(X,istart,(i+istart-1)) |
|
395 ! r = isize |
|
396 r = 1 |
|
397 Call SwapDown(X,r,i-1,istart,keyfield) |
|
398 ! |
|
399 ! print *, "P after Swapdown is " |
|
400 ! do j=istart,isize |
|
401 ! print *, j, X(j,keyfield) |
|
402 ! end do |
|
403 end do |
|
404 ! |
|
405 return |
|
406 end subroutine HeapSort |
|
407 ! |
|
408 !------------------------------------------------------------------------------ |
|
409 ! |
|
410 subroutine XPMMap(N,II) |
|
411 ! |
|
412 ! output routine for xpm X-window display tool |
|
413 ! |
|
414 ! |
|
415 implicit none |
|
416 integer(i4), intent(in) :: N |
|
417 integer(i4), intent(in) :: II(0:maxN,0:maxN) |
|
418 ! local variables |
|
419 integer(i4) :: itemp(0:maxN) |
|
420 integer(i4) :: i, j ! loop indices |
|
421 character (len=40) :: fmt |
|
422 ! |
|
423 ! |
|
424 ! output map |
|
425 print *, "writing final map to xpm format file landscape.xpm" |
|
426 ! |
|
427 open (20,file='landscape.xpm',status='unknown') |
|
428 ! |
|
429 write (20, 100) |
|
430 100 format ('/* XPM */',/'static char *realizer[]', & |
|
431 ' = {',/'/* cols rows ncolors chars_per_pixel', & |
|
432 ' */') |
|
433 write (20, 102) N+1, N+1 |
|
434 102 format ('"',i6,1x,i6,1x,'11 1 0 0",',/'/* colors */') |
|
435 write (20, 104) |
|
436 104 format ('"0 c white m white', & |
|
437 '",',/'"1 s Decid c BurlyWood m black', & |
|
438 '",',/'"2 s Mixed c DarkGoldenrod3 m gray40', & |
|
439 '",',/'"3 s Everg c MistyRose3 m gray50', & |
|
440 '",',/'"4 s Trans c LightGoldenrod4 m gray60', & |
|
441 '",',/'"5 s Urban c DarkGoldenrod1 m gray70', & |
|
442 '",',/'"6 s Water c gold m gray80', & |
|
443 '",',/'"7 s Maint c LightGoldenrod3 m white', & |
|
444 '",',/'"8 s Cropl c SaddleBrown m lightgray', & |
|
445 '",',/'"9 s Lawng c RosyBrown m gray', & |
|
446 '",',/'"* c chocolate m black",',/'/* pixels */') |
|
447 ! |
|
448 write(fmt,'("(",I0,"(i1)$)")') N+1 |
|
449 do i= 0, N |
|
450 do j = 0, N |
|
451 ! itemp(j) = II(i,j) |
|
452 itemp(j) = II(j,i) |
|
453 end do |
|
454 write (20,'(1A$)') '"' |
|
455 ! write (20,108) (itemp(j), j=0,N) |
|
456 write (20,fmt) (itemp(j), j=0,N) |
|
457 write (20,'(2A)') '",' |
|
458 !108 format (' "',<N+1>(i1),'",') |
|
459 end do |
|
460 ! |
|
461 ! |
|
462 write (20,112) |
|
463 112 format ('};') |
|
464 ! |
|
465 return |
|
466 end subroutine XPMMap |
|
467 ! |
|
468 !------------------------------------------------------------------------------ |
|
469 ! |
|
470 subroutine InputMap(mapfile,P,N,k,invert) |
|
471 ! |
|
472 ! input map from a data file |
|
473 ! |
|
474 implicit none |
|
475 character (len=60), intent(in) :: mapfile |
|
476 integer(i4), intent(inout) :: P(maxN*maxN,0:maxsp+3) |
|
477 integer(i4), intent(in) :: N |
|
478 integer(i4), intent(in) :: k |
|
479 logical, intent(in) :: invert |
|
480 ! local variables |
|
481 integer(i4) :: i, j, iv ! loop indices |
|
482 real(r4) :: X(0:maxN,0:maxN) |
|
483 real(r4) :: xmin |
|
484 real(r4) :: xmax |
|
485 real(r4) :: xrange |
|
486 real(r4) :: temp |
|
487 ! |
|
488 open (14, file=mapfile, status='old') |
|
489 ! |
|
490 xmin = 1e30 |
|
491 xmax = -1e30 |
|
492 do j = 0, N |
|
493 read (14,*) ( X(i,j), i = 0, N) |
|
494 ! read (14,118) X(i,j) |
|
495 !118 format (<(N+1)*(N+1)>i5, 1x) |
|
496 do i = 0, N |
|
497 if (X(i,j) .lt. xmin) xmin = X(i,j) |
|
498 if (X(i,j) .gt. xmax) xmax = X(i,j) |
|
499 end do |
|
500 end do |
|
501 ! |
|
502 ! invert if necessary and scale between 0 and 32767 |
|
503 xrange = xmax - xmin |
|
504 do j = 0, N |
|
505 do i = 0, N |
|
506 if (invert) then |
|
507 temp = ((-X(i,j) + xmax) / xrange) * 32767 |
|
508 else |
|
509 temp = ((X(i,j) - xmin) / xrange) * 32767 |
|
510 end if |
|
511 ! |
|
512 iv = (i*(N+1))+j+1 |
|
513 P(iv,k) = nint(temp) |
|
514 P(iv,ihabs+1) = i |
|
515 P(iv,ihabs+2) = j |
|
516 end do |
|
517 end do |
|
518 ! |
|
519 close(14) |
|
520 rewind(14) |
|
521 ! |
|
522 return |
|
523 end subroutine InputMap |
|
524 ! |
|
525 !------------------------------------------------------------------------------ |
|
526 ! |
|
527 subroutine SD2(Heap,r,isize,istart) |
|
528 implicit none |
|
529 integer(i4), intent(inout) :: Heap(maxN*maxN) |
|
530 integer(i4), intent(inout) :: r |
|
531 integer(i4), intent(in) :: isize |
|
532 integer(i4), intent(in) :: istart |
|
533 ! local variables |
|
534 integer(i4) :: ipos |
|
535 integer(i4) :: child |
|
536 logical :: done |
|
537 ! |
|
538 ! print *, 'entering SwapDown' |
|
539 ! print *, 'params are ', r, isize, istart |
|
540 ipos = istart - 1 |
|
541 done = .false. |
|
542 child = 2 * r |
|
543 do while (.not. done .and. child .le. isize) |
|
544 ! find the largest child |
|
545 if (child .lt. isize .and. Heap(child+ipos) .lt. & |
|
546 Heap(child+ipos+1)) child = child + 1 |
|
547 ! interchange node and largest child if necessary |
|
548 ! and move down to the next subtree |
|
549 if (Heap(r+ipos) .lt. Heap(child+ipos)) then |
|
550 Call I2(Heap,r+ipos,child+ipos) |
|
551 r = child |
|
552 child = 2 * child |
|
553 else |
|
554 done = .true. |
|
555 end if |
|
556 ! |
|
557 end do |
|
558 ! |
|
559 return |
|
560 end subroutine SD2 |
|
561 ! |
|
562 !------------------------------------------------------------------------------ |
|
563 ! |
|
564 subroutine I2(X,m,n) |
|
565 implicit none |
|
566 integer(i4), intent(inout) :: X(maxN*maxN) |
|
567 integer(i4), intent(in) :: m |
|
568 integer(i4), intent(in) :: n |
|
569 ! local variables |
|
570 integer(i4) :: temp |
|
571 ! |
|
572 temp = X(m) |
|
573 X(m) = X(n) |
|
574 X(n) = temp |
|
575 ! |
|
576 return |
|
577 end subroutine I2 |
|
578 ! |
|
579 !------------------------------------------------------------------------------ |
|
580 ! |
|
581 subroutine HS2(X,isize,istart) |
|
582 implicit none |
|
583 integer(i4), intent(inout) :: X(maxN*maxN) |
|
584 integer(i4), intent(in) :: isize |
|
585 integer(i4), intent(in) :: istart |
|
586 ! local variables |
|
587 integer(i4) :: i, j ! loop indices |
|
588 integer(i4) :: r |
|
589 integer(i4) :: tempr |
|
590 ! |
|
591 ! print *, "isize in HeapSort is ", isize |
|
592 ! print *, "istart in HeapSort is ", istart |
|
593 ! heapify |
|
594 do r = isize/2, 1, -1 |
|
595 tempr = r |
|
596 Call SD2(X,tempr,isize,istart) |
|
597 !Call SD2(X,r,isize,istart) |
|
598 !r = tempr |
|
599 end do |
|
600 ! |
|
601 ! print*, "P after heapify is " |
|
602 ! do i=istart,isize |
|
603 ! print *, i, X(i) |
|
604 ! end do |
|
605 ! |
|
606 ! repeatedly put largest item in root at end of list, |
|
607 ! prune it from the tree |
|
608 ! and apply SwapDown to the rest of the tree |
|
609 ! |
|
610 do i = isize, 2, -1 |
|
611 Call I2(X,istart,(i+istart-1)) |
|
612 ! r = isize |
|
613 r = 1 |
|
614 Call SD2(X,r,i-1,istart) |
|
615 ! |
|
616 ! print*, "P after Swapdown is " |
|
617 ! do j=istart,isize |
|
618 ! print *, j, X(j) |
|
619 ! end do |
|
620 |
|
621 end do |
|
622 ! |
|
623 return |
|
624 end subroutine HS2 |
|
625 ! |
|
626 !------------------------------------------------------------------------------ |
|
627 ! |
|
628 subroutine XPMTies(P) |
|
629 ! |
|
630 ! output routine for xpm X-window display tool |
|
631 ! |
|
632 ! write xpm format map with grayshades |
|
633 ! |
|
634 implicit none |
|
635 integer(i4), intent(in) :: P(maxN*maxN,0:maxsp+3) |
|
636 ! local variables |
|
637 integer(i4) :: itemp(0:maxN,0:maxN) |
|
638 integer(i4) :: i, j, m ! loop indices |
|
639 real(r4) :: scalelog |
|
640 character (len=40) :: fmt |
|
641 ! |
|
642 open (40,file='ties.xpm',status='unknown') |
|
643 ! |
|
644 print 99 |
|
645 99 format ('XPM map file ties.xpm being written') |
|
646 write (40, 100) |
|
647 100 format ('/* XPM */',/'static char *tiesmap[]', & |
|
648 ' = {',/'/* cols rows ncolors chars_per_pixel', & |
|
649 ' */') |
|
650 write (40, 102) irowcols+1, irowcols+1 |
|
651 102 format (' "',i6,1x,i6,1x,'108 3 0 0",',/'/* colors */') |
|
652 write (40, 104) |
|
653 104 format ('" 0 s 0.0prob c white g4 white g white m white', & |
|
654 '",',/'" 1 s 0.01prob c gray99 g4 gray99 g gray99 m gray99', & |
|
655 '",',/'" 2 s 0.02prob c gray98 g4 gray98 g gray98 m gray98', & |
|
656 '",',/'" 3 s 0.03prob c gray97 g4 gray97 g gray97 m gray97', & |
|
657 '",',/'" 4 s 0.04prob c gray96 g4 gray96 g gray96 m gray96', & |
|
658 '",',/'" 5 s 0.05prob c gray95 g4 gray95 g gray95 m gray95', & |
|
659 '",',/'" 6 s 0.06prob c gray94 g4 gray94 g gray94 m gray94', & |
|
660 '",',/'" 7 s 0.07prob c gray93 g4 gray93 g gray93 m gray93', & |
|
661 '",',/'" 8 s 0.08prob c gray92 g4 gray92 g gray92 m gray92', & |
|
662 '",',/'" 9 s 0.09prob c gray91 g4 gray91 g gray91 m gray91', & |
|
663 '",',/'" 10 s 0.10prob c gray90 g4 gray90 g gray90 m gray90', & |
|
664 '",',/'" 11 s 0.11prob c gray89 g4 gray89 g gray89 m gray89', & |
|
665 '",',/'" 12 s 0.12prob c gray88 g4 gray88 g gray88 m gray88', & |
|
666 '",',/'" 13 s 0.13prob c gray87 g4 gray87 g gray87 m gray87', & |
|
667 '",',/'" 14 s 0.14prob c gray86 g4 gray86 g gray86 m gray86', & |
|
668 '",',/'" 15 s 0.15prob c gray85 g4 gray85 g gray85 m gray85', & |
|
669 '",',/'" 16 s 0.16prob c gray84 g4 gray84 g gray84 m gray84', & |
|
670 '",',/'" 17 s 0.17prob c gray83 g4 gray83 g gray83 m gray83', & |
|
671 '",',/'" 18 s 0.18prob c gray82 g4 gray82 g gray82 m gray82', & |
|
672 '",',/'" 19 s 0.19prob c gray81 g4 gray81 g gray81 m gray81', & |
|
673 '",',/'" 20 s 0.20prob c gray80 g4 gray80 g gray80 m gray80', & |
|
674 '",',/'" 21 s 0.21prob c gray79 g4 gray79 g gray79 m gray79', & |
|
675 '",',/'" 22 s 0.22prob c gray78 g4 gray78 g gray78 m gray78', & |
|
676 '",',/'" 23 s 0.23prob c gray77 g4 gray77 g gray77 m gray77', & |
|
677 '",',/'" 24 s 0.24prob c gray76 g4 gray76 g gray76 m gray76', & |
|
678 '",',/'" 25 s 0.25prob c gray75 g4 gray75 g gray75 m gray75', & |
|
679 '",',/'" 26 s 0.26prob c gray74 g4 gray74 g gray74 m gray74', & |
|
680 '",',/'" 27 s 0.27prob c gray73 g4 gray73 g gray73 m gray73', & |
|
681 '",',/'" 28 s 0.28prob c gray72 g4 gray72 g gray72 m gray72', & |
|
682 '",',/'" 29 s 0.29prob c gray71 g4 gray71 g gray71 m gray71', & |
|
683 '",',/'" 30 s 0.30prob c gray70 g4 gray70 g gray70 m gray70', & |
|
684 '",',/'" 31 s 0.31prob c gray69 g4 gray69 g gray69 m gray69', & |
|
685 '",',/'" 32 s 0.32prob c gray68 g4 gray68 g gray68 m gray68', & |
|
686 '",',/'" 33 s 0.33prob c gray67 g4 gray67 g gray67 m gray67', & |
|
687 '",',/'" 34 s 0.34prob c gray66 g4 gray66 g gray66 m gray66', & |
|
688 '",',/'" 35 s 0.35prob c gray65 g4 gray65 g gray65 m gray65', & |
|
689 '",',/'" 36 s 0.36prob c gray64 g4 gray64 g gray64 m gray64', & |
|
690 '",',/'" 37 s 0.37prob c gray63 g4 gray63 g gray63 m gray63', & |
|
691 '",',/'" 38 s 0.38prob c gray62 g4 gray62 g gray62 m gray62', & |
|
692 '",',/'" 39 s 0.39prob c gray61 g4 gray61 g gray61 m gray61', & |
|
693 '",',/'" 40 s 0.40prob c gray60 g4 gray60 g gray60 m gray60', & |
|
694 '",',/'" 41 s 0.41prob c gray59 g4 gray59 g gray59 m gray59', & |
|
695 '",',/'" 42 s 0.42prob c gray58 g4 gray58 g gray58 m gray58', & |
|
696 '",',/'" 43 s 0.43prob c gray57 g4 gray57 g gray57 m gray57', & |
|
697 '",',/'" 44 s 0.44prob c gray56 g4 gray56 g gray56 m gray56', & |
|
698 '",',/'" 45 s 0.45prob c gray55 g4 gray55 g gray55 m gray55', & |
|
699 '",',/'" 46 s 0.46prob c gray54 g4 gray54 g gray54 m gray54', & |
|
700 '",',/'" 47 s 0.47prob c gray53 g4 gray53 g gray53 m gray53', & |
|
701 '",',/'" 48 s 0.48prob c gray52 g4 gray52 g gray52 m gray52', & |
|
702 '",',/'" 49 s 0.49prob c gray51 g4 gray51 g gray51 m gray51', & |
|
703 '",',/'" 50 s 0.50prob c gray50 g4 gray50 g gray50 m gray50', & |
|
704 '",',/'" 51 s 0.51prob c gray49 g4 gray49 g gray49 m gray49', & |
|
705 '",',/'" 52 s 0.52prob c gray48 g4 gray48 g gray48 m gray48', & |
|
706 '",',/'" 53 s 0.53prob c gray47 g4 gray47 g gray47 m gray47', & |
|
707 '",',/'" 54 s 0.54prob c gray46 g4 gray46 g gray46 m gray46', & |
|
708 '",',/'" 55 s 0.55prob c gray45 g4 gray45 g gray45 m gray45', & |
|
709 '",',/'" 56 s 0.56prob c gray44 g4 gray44 g gray44 m gray44', & |
|
710 '",',/'" 57 s 0.57prob c gray43 g4 gray43 g gray43 m gray43', & |
|
711 '",',/'" 58 s 0.58prob c gray42 g4 gray42 g gray42 m gray42', & |
|
712 '",',/'" 59 s 0.59prob c gray41 g4 gray41 g gray41 m gray41', & |
|
713 '",',/'" 60 s 0.60prob c gray40 g4 gray40 g gray40 m gray40', & |
|
714 '",',/'" 61 s 0.61prob c gray39 g4 gray39 g gray39 m gray39', & |
|
715 '",',/'" 62 s 0.62prob c gray38 g4 gray38 g gray38 m gray38', & |
|
716 '",',/'" 63 s 0.63prob c gray37 g4 gray37 g gray37 m gray37', & |
|
717 '",',/'" 64 s 0.64prob c gray36 g4 gray36 g gray36 m gray36', & |
|
718 '",',/'" 65 s 0.65prob c gray35 g4 gray35 g gray35 m gray35', & |
|
719 '",',/'" 66 s 0.66prob c gray34 g4 gray34 g gray34 m gray34', & |
|
720 '",',/'" 67 s 0.67prob c gray33 g4 gray33 g gray33 m gray33', & |
|
721 '",',/'" 68 s 0.68prob c gray32 g4 gray32 g gray32 m gray32', & |
|
722 '",',/'" 69 s 0.69prob c gray31 g4 gray31 g gray31 m gray31', & |
|
723 '",',/'" 70 s 0.70prob c gray30 g4 gray30 g gray30 m gray30', & |
|
724 '",',/'" 71 s 0.71prob c gray29 g4 gray29 g gray29 m gray29', & |
|
725 '",',/'" 72 s 0.72prob c gray28 g4 gray28 g gray28 m gray28', & |
|
726 '",',/'" 73 s 0.73prob c gray27 g4 gray27 g gray27 m gray27', & |
|
727 '",',/'" 74 s 0.74prob c gray26 g4 gray26 g gray26 m gray26', & |
|
728 '",',/'" 75 s 0.75prob c gray25 g4 gray25 g gray25 m gray25', & |
|
729 '",',/'" 76 s 0.76prob c gray24 g4 gray24 g gray24 m gray24', & |
|
730 '",',/'" 77 s 0.77prob c gray23 g4 gray23 g gray23 m gray23', & |
|
731 '",',/'" 78 s 0.78prob c gray22 g4 gray22 g gray22 m gray22', & |
|
732 '",',/'" 79 s 0.79prob c gray21 g4 gray21 g gray21 m gray21', & |
|
733 '",',/'" 80 s 0.80prob c gray20 g4 gray20 g gray20 m gray20', & |
|
734 '",',/'" 81 s 0.81prob c gray19 g4 gray19 g gray19 m gray19', & |
|
735 '",',/'" 82 s 0.82prob c gray18 g4 gray18 g gray18 m gray18', & |
|
736 '",',/'" 83 s 0.83prob c gray17 g4 gray17 g gray17 m gray17', & |
|
737 '",',/'" 84 s 0.84prob c gray16 g4 gray16 g gray16 m gray16', & |
|
738 '",',/'" 85 s 0.85prob c gray15 g4 gray15 g gray15 m gray15', & |
|
739 '",',/'" 86 s 0.86prob c gray14 g4 gray14 g gray14 m gray14', & |
|
740 '",',/'" 87 s 0.87prob c gray13 g4 gray13 g gray13 m gray13', & |
|
741 '",',/'" 88 s 0.88prob c gray12 g4 gray12 g gray12 m gray12', & |
|
742 '",',/'" 89 s 0.89prob c gray11 g4 gray11 g gray11 m gray11', & |
|
743 '",',/'" 90 s 0.90prob c gray10 g4 gray10 g gray10 m gray10', & |
|
744 '",',/'" 91 s 0.91prob c gray9 g4 gray9 g gray9 m gray9', & |
|
745 '",',/'" 92 s 0.92prob c gray8 g4 gray8 g gray8 m gray8', & |
|
746 '",',/'" 93 s 0.93prob c gray7 g4 gray7 g gray7 m gray7', & |
|
747 '",',/'" 94 s 0.94prob c gray6 g4 gray6 g gray6 m gray6', & |
|
748 '",',/'" 95 s 0.95prob c gray5 g4 gray5 g gray5 m gray5', & |
|
749 '",',/'" 96 s 0.96prob c gray4 g4 gray4 g gray4 m gray4', & |
|
750 '",',/'" 97 s 0.97prob c gray3 g4 gray3 g gray3 m gray3', & |
|
751 '",',/'" 98 s 0.98prob c gray2 g4 gray2 g gray2 m gray2', & |
|
752 '",',/'" 99 s 0.99prob c gray1 g4 gray1 g gray1 m gray1', & |
|
753 '",',/'"100 s 1.00prob c black g4 black g black m black', & |
|
754 '",',/'"101 c white m white s NotFlammable', & |
|
755 '",',/'"102 c GreenYellow m gray70 s LP0', & |
|
756 '",',/'"103 c LawnGreen m gray60 s LP1', & |
|
757 '",',/'"104 c LimeGreen m gray50 s LP2', & |
|
758 '",',/'"105 c ForestGreen m gray40 s LP3', & |
|
759 '",',/'"106 c yellow m gray80 s NF', & |
|
760 '",',/'"107 c blue m black s water ",', & |
|
761 /'/* pixels */') |
|
762 ! |
|
763 ! calculate a log scale for gray levels |
|
764 ! use this scale log for uniform scaling across realizations |
|
765 ! scalelog = 10**(log10(100.)/32767.) |
|
766 ! scalelog for simple summed priority |
|
767 ! scale max gray to ihabs-isupply-1 full-on categories |
|
768 ! (= ihabs-isupply-1 x 32767) |
|
769 !c scalelog = 10**(log10(100.)/(32767. * (ihabs-isupply-1))) |
|
770 ! use this scalelog for maximum priority resolution w/in one map |
|
771 scalelog = 10**(log10(100.) / real(maxgray,kind(scalelog))) |
|
772 ! print *, scalelog |
|
773 ! |
|
774 ! start at bottom of list (highest priority ties) |
|
775 do m = numcells, itiedsplit+1, -1 ! over all tied |
|
776 ! itemp(P(m,ihabs+1),P(m,ihabs+2)) = |
|
777 ! & (real(P(m,ihabs+3),kind(itemp))/32767)*100 |
|
778 ! |
|
779 ! print *, scalelog**P(m,ihabs+3) |
|
780 ! print *, P(m,ihabs+1),P(m,ihabs+2),P(m,ihabs+3) |
|
781 itemp(P(m,ihabs+1),P(m,ihabs+2)) = scalelog**P(m,ihabs+3) |
|
782 ! itemp(P(m,ihabs+1),P(m,ihabs+2)) = 100 |
|
783 ! all priorities greater than color index 100 set black |
|
784 if (itemp(P(m,ihabs+1),P(m,ihabs+2)) .gt. 100) & |
|
785 itemp(P(m,ihabs+1),P(m,ihabs+2)) = 100 |
|
786 end do |
|
787 ! |
|
788 ! color used blanks yellow |
|
789 do m = mtsplit,mtresolved+1,-1 |
|
790 itemp(P(m,ihabs+1),P(m,ihabs+2)) = 106 |
|
791 end do |
|
792 ! |
|
793 write (fmt,'("(",I0,"(i3)$)")') irowcols+1 |
|
794 do i = 0, irowcols |
|
795 write (40,'(1A$)') '"' |
|
796 !write (40,108) (itemp(j,i), j=0, irowcols) |
|
797 write (40,fmt) (itemp(j,i), j=0, irowcols) |
|
798 write (40,'(2A)') '",' |
|
799 !108 format ('"',<irowcols+1>(i3),'",') |
|
800 end do |
|
801 ! |
|
802 write (40,112) |
|
803 112 format ('};') |
|
804 close (40) |
|
805 ! |
|
806 return |
|
807 end subroutine XPMTies |
|
808 ! |
|
809 !------------------------------------------------------------------------------ |
|
810 ! |
|
811 subroutine XPMRelief(P,iz,iflood) |
|
812 ! |
|
813 ! output routine for xpm X-window display tool |
|
814 ! |
|
815 ! write xpm format aspect map of painted fractal |
|
816 ! spatial probability surfaces with grayshades |
|
817 ! |
|
818 implicit none |
|
819 integer(i4), intent(in) :: P(maxN*maxN,0:maxsp+3) |
|
820 integer(i4), intent(in) :: iz |
|
821 integer(i4), intent(in) :: iflood |
|
822 ! local variables |
|
823 integer(i4) :: itemp(0:maxN,0:maxN) |
|
824 integer(i4) :: idispvec(0:maxN) |
|
825 integer(i4) :: i, j, k, l, m ! loop indices |
|
826 real(r4) :: alt ! altitude |
|
827 real(r4) :: az ! azimuth |
|
828 real(r4) :: x, y |
|
829 real(r4) :: slope ! slope |
|
830 real(r4) :: aspect ! aspect |
|
831 real(r4) :: c |
|
832 character (len=15) :: outfile(0:maxsp+1) |
|
833 character (len=40) :: fmt |
|
834 ! |
|
835 ! set filenames for the first 10 categories |
|
836 outfile(0) = 'probmap0.xpm ' |
|
837 outfile(1) = 'probmap1.xpm ' |
|
838 outfile(2) = 'probmap2.xpm ' |
|
839 outfile(3) = 'probmap3.xpm ' |
|
840 outfile(4) = 'probmap4.xpm ' |
|
841 outfile(5) = 'probmap5.xpm ' |
|
842 outfile(6) = 'probmap6.xpm ' |
|
843 outfile(7) = 'probmap7.xpm ' |
|
844 outfile(8) = 'probmap8.xpm ' |
|
845 outfile(9) = 'probmap9.xpm ' |
|
846 outfile(10) = 'probmap10.xpm ' |
|
847 ! |
|
848 ! |
|
849 print 99, outfile(iz) |
|
850 99 format ('XPM map file ', a15 'being written') |
|
851 ! |
|
852 ! open the right file for xpm output |
|
853 open (50,file=outfile(iz),status='unknown') |
|
854 ! |
|
855 write (50, 100) |
|
856 100 format ('/* XPM */',/'static char *fracasp[]', & |
|
857 ' = {',/'/* cols rows ncolors chars_per_pixel', ' */') |
|
858 write (50, 102) irowcols+1, irowcols+1 |
|
859 102 format ('"',i6,1x,i6,1x,'108 3 0 0",',/'/* colors */') |
|
860 write (50, 104) |
|
861 104 format ('" 0 s 0.0prob c white g4 white g white m white', & |
|
862 '",',/'" 1 s 0.01prob c gray99 g4 gray99 g gray99 m gray99', & |
|
863 '",',/'" 2 s 0.02prob c gray98 g4 gray98 g gray98 m gray98', & |
|
864 '",',/'" 3 s 0.03prob c gray97 g4 gray97 g gray97 m gray97', & |
|
865 '",',/'" 4 s 0.04prob c gray96 g4 gray96 g gray96 m gray96', & |
|
866 '",',/'" 5 s 0.05prob c gray95 g4 gray95 g gray95 m gray95', & |
|
867 '",',/'" 6 s 0.06prob c gray94 g4 gray94 g gray94 m gray94', & |
|
868 '",',/'" 7 s 0.07prob c gray93 g4 gray93 g gray93 m gray93', & |
|
869 '",',/'" 8 s 0.08prob c gray92 g4 gray92 g gray92 m gray92', & |
|
870 '",',/'" 9 s 0.09prob c gray91 g4 gray91 g gray91 m gray91', & |
|
871 '",',/'" 10 s 0.10prob c gray90 g4 gray90 g gray90 m gray90', & |
|
872 '",',/'" 11 s 0.11prob c gray89 g4 gray89 g gray89 m gray89', & |
|
873 '",',/'" 12 s 0.12prob c gray88 g4 gray88 g gray88 m gray88', & |
|
874 '",',/'" 13 s 0.13prob c gray87 g4 gray87 g gray87 m gray87', & |
|
875 '",',/'" 14 s 0.14prob c gray86 g4 gray86 g gray86 m gray86', & |
|
876 '",',/'" 15 s 0.15prob c gray85 g4 gray85 g gray85 m gray85', & |
|
877 '",',/'" 16 s 0.16prob c gray84 g4 gray84 g gray84 m gray84', & |
|
878 '",',/'" 17 s 0.17prob c gray83 g4 gray83 g gray83 m gray83', & |
|
879 '",',/'" 18 s 0.18prob c gray82 g4 gray82 g gray82 m gray82', & |
|
880 '",',/'" 19 s 0.19prob c gray81 g4 gray81 g gray81 m gray81', & |
|
881 '",',/'" 20 s 0.20prob c gray80 g4 gray80 g gray80 m gray80', & |
|
882 '",',/'" 21 s 0.21prob c gray79 g4 gray79 g gray79 m gray79', & |
|
883 '",',/'" 22 s 0.22prob c gray78 g4 gray78 g gray78 m gray78', & |
|
884 '",',/'" 23 s 0.23prob c gray77 g4 gray77 g gray77 m gray77', & |
|
885 '",',/'" 24 s 0.24prob c gray76 g4 gray76 g gray76 m gray76', & |
|
886 '",',/'" 25 s 0.25prob c gray75 g4 gray75 g gray75 m gray75', & |
|
887 '",',/'" 26 s 0.26prob c gray74 g4 gray74 g gray74 m gray74', & |
|
888 '",',/'" 27 s 0.27prob c gray73 g4 gray73 g gray73 m gray73', & |
|
889 '",',/'" 28 s 0.28prob c gray72 g4 gray72 g gray72 m gray72', & |
|
890 '",',/'" 29 s 0.29prob c gray71 g4 gray71 g gray71 m gray71', & |
|
891 '",',/'" 30 s 0.30prob c gray70 g4 gray70 g gray70 m gray70', & |
|
892 '",',/'" 31 s 0.31prob c gray69 g4 gray69 g gray69 m gray69', & |
|
893 '",',/'" 32 s 0.32prob c gray68 g4 gray68 g gray68 m gray68', & |
|
894 '",',/'" 33 s 0.33prob c gray67 g4 gray67 g gray67 m gray67', & |
|
895 '",',/'" 34 s 0.34prob c gray66 g4 gray66 g gray66 m gray66', & |
|
896 '",',/'" 35 s 0.35prob c gray65 g4 gray65 g gray65 m gray65', & |
|
897 '",',/'" 36 s 0.36prob c gray64 g4 gray64 g gray64 m gray64', & |
|
898 '",',/'" 37 s 0.37prob c gray63 g4 gray63 g gray63 m gray63', & |
|
899 '",',/'" 38 s 0.38prob c gray62 g4 gray62 g gray62 m gray62', & |
|
900 '",',/'" 39 s 0.39prob c gray61 g4 gray61 g gray61 m gray61', & |
|
901 '",',/'" 40 s 0.40prob c gray60 g4 gray60 g gray60 m gray60', & |
|
902 '",',/'" 41 s 0.41prob c gray59 g4 gray59 g gray59 m gray59', & |
|
903 '",',/'" 42 s 0.42prob c gray58 g4 gray58 g gray58 m gray58', & |
|
904 '",',/'" 43 s 0.43prob c gray57 g4 gray57 g gray57 m gray57', & |
|
905 '",',/'" 44 s 0.44prob c gray56 g4 gray56 g gray56 m gray56', & |
|
906 '",',/'" 45 s 0.45prob c gray55 g4 gray55 g gray55 m gray55', & |
|
907 '",',/'" 46 s 0.46prob c gray54 g4 gray54 g gray54 m gray54', & |
|
908 '",',/'" 47 s 0.47prob c gray53 g4 gray53 g gray53 m gray53', & |
|
909 '",',/'" 48 s 0.48prob c gray52 g4 gray52 g gray52 m gray52', & |
|
910 '",',/'" 49 s 0.49prob c gray51 g4 gray51 g gray51 m gray51', & |
|
911 '",',/'" 50 s 0.50prob c gray50 g4 gray50 g gray50 m gray50', & |
|
912 '",',/'" 51 s 0.51prob c gray49 g4 gray49 g gray49 m gray49', & |
|
913 '",',/'" 52 s 0.52prob c gray48 g4 gray48 g gray48 m gray48', & |
|
914 '",',/'" 53 s 0.53prob c gray47 g4 gray47 g gray47 m gray47', & |
|
915 '",',/'" 54 s 0.54prob c gray46 g4 gray46 g gray46 m gray46', & |
|
916 '",',/'" 55 s 0.55prob c gray45 g4 gray45 g gray45 m gray45', & |
|
917 '",',/'" 56 s 0.56prob c gray44 g4 gray44 g gray44 m gray44', & |
|
918 '",',/'" 57 s 0.57prob c gray43 g4 gray43 g gray43 m gray43', & |
|
919 '",',/'" 58 s 0.58prob c gray42 g4 gray42 g gray42 m gray42', & |
|
920 '",',/'" 59 s 0.59prob c gray41 g4 gray41 g gray41 m gray41', & |
|
921 '",',/'" 60 s 0.60prob c gray40 g4 gray40 g gray40 m gray40', & |
|
922 '",',/'" 61 s 0.61prob c gray39 g4 gray39 g gray39 m gray39', & |
|
923 '",',/'" 62 s 0.62prob c gray38 g4 gray38 g gray38 m gray38', & |
|
924 '",',/'" 63 s 0.63prob c gray37 g4 gray37 g gray37 m gray37', & |
|
925 '",',/'" 64 s 0.64prob c gray36 g4 gray36 g gray36 m gray36', & |
|
926 '",',/'" 65 s 0.65prob c gray35 g4 gray35 g gray35 m gray35', & |
|
927 '",',/'" 66 s 0.66prob c gray34 g4 gray34 g gray34 m gray34', & |
|
928 '",',/'" 67 s 0.67prob c gray33 g4 gray33 g gray33 m gray33', & |
|
929 '",',/'" 68 s 0.68prob c gray32 g4 gray32 g gray32 m gray32', & |
|
930 '",',/'" 69 s 0.69prob c gray31 g4 gray31 g gray31 m gray31', & |
|
931 '",',/'" 70 s 0.70prob c gray30 g4 gray30 g gray30 m gray30', & |
|
932 '",',/'" 71 s 0.71prob c gray29 g4 gray29 g gray29 m gray29', & |
|
933 '",',/'" 72 s 0.72prob c gray28 g4 gray28 g gray28 m gray28', & |
|
934 '",',/'" 73 s 0.73prob c gray27 g4 gray27 g gray27 m gray27', & |
|
935 '",',/'" 74 s 0.74prob c gray26 g4 gray26 g gray26 m gray26', & |
|
936 '",',/'" 75 s 0.75prob c gray25 g4 gray25 g gray25 m gray25', & |
|
937 '",',/'" 76 s 0.76prob c gray24 g4 gray24 g gray24 m gray24', & |
|
938 '",',/'" 77 s 0.77prob c gray23 g4 gray23 g gray23 m gray23', & |
|
939 '",',/'" 78 s 0.78prob c gray22 g4 gray22 g gray22 m gray22', & |
|
940 '",',/'" 79 s 0.79prob c gray21 g4 gray21 g gray21 m gray21', & |
|
941 '",',/'" 80 s 0.80prob c gray20 g4 gray20 g gray20 m gray20', & |
|
942 '",',/'" 81 s 0.81prob c gray19 g4 gray19 g gray19 m gray19', & |
|
943 '",',/'" 82 s 0.82prob c gray18 g4 gray18 g gray18 m gray18', & |
|
944 '",',/'" 83 s 0.83prob c gray17 g4 gray17 g gray17 m gray17', & |
|
945 '",',/'" 84 s 0.84prob c gray16 g4 gray16 g gray16 m gray16', & |
|
946 '",',/'" 85 s 0.85prob c gray15 g4 gray15 g gray15 m gray15', & |
|
947 '",',/'" 86 s 0.86prob c gray14 g4 gray14 g gray14 m gray14', & |
|
948 '",',/'" 87 s 0.87prob c gray13 g4 gray13 g gray13 m gray13', & |
|
949 '",',/'" 88 s 0.88prob c gray12 g4 gray12 g gray12 m gray12', & |
|
950 '",',/'" 89 s 0.89prob c gray11 g4 gray11 g gray11 m gray11', & |
|
951 '",',/'" 90 s 0.90prob c gray10 g4 gray10 g gray10 m gray10', & |
|
952 '",',/'" 91 s 0.91prob c gray9 g4 gray9 g gray9 m gray9', & |
|
953 '",',/'" 92 s 0.92prob c gray8 g4 gray8 g gray8 m gray8', & |
|
954 '",',/'" 93 s 0.93prob c gray7 g4 gray7 g gray7 m gray7', & |
|
955 '",',/'" 94 s 0.94prob c gray6 g4 gray6 g gray6 m gray6', & |
|
956 '",',/'" 95 s 0.95prob c gray5 g4 gray5 g gray5 m gray5', & |
|
957 '",',/'" 96 s 0.96prob c gray4 g4 gray4 g gray4 m gray4', & |
|
958 '",',/'" 97 s 0.97prob c gray3 g4 gray3 g gray3 m gray3', & |
|
959 '",',/'" 98 s 0.98prob c gray2 g4 gray2 g gray2 m gray2', & |
|
960 '",',/'" 99 s 0.99prob c gray1 g4 gray1 g gray1 m gray1', & |
|
961 '",',/'"100 s 1.00prob c black g4 black g black m black', & |
|
962 '",',/'"101 c white m white s NotFlammable', & |
|
963 '",',/'"102 c GreenYellow m gray70 s LP0', & |
|
964 '",',/'"103 c LawnGreen m gray60 s LP1', & |
|
965 '",',/'"104 c LimeGreen m gray50 s LP2', & |
|
966 '",',/'"105 c ForestGreen m gray40 s LP3', & |
|
967 '",',/'"106 c yellow m gray80 s NF', & |
|
968 '",',/'"107 c blue m black s water ",', & |
|
969 /'/* pixels */') |
|
970 ! |
|
971 ! |
|
972 ! load fractal probability terrain into itemp array |
|
973 do m = 1, numcells ! over all cells |
|
974 itemp(P(m,ihabs+1),P(m,ihabs+2)) = P(m,iz) |
|
975 end do |
|
976 ! |
|
977 ! |
|
978 ! equations for slope and aspect are from: |
|
979 ! Horn, B.K.P., 1981. Hill Shading and the Reflectance Map. |
|
980 ! Proceedings of the I.E.E.E., 69(1):14-47. |
|
981 ! |
|
982 ! calculate and display aspect map |
|
983 ! set sun altitude above horizon in degrees (0 to 90) |
|
984 alt = 40. |
|
985 ! set sun azimuth in degrees west of north (-1 to 360) |
|
986 az = 105. |
|
987 k = 0 |
|
988 do j = 0,irowcols |
|
989 do i = 0,irowcols |
|
990 ! |
|
991 x = (itemp(i-1,j-1)+2.*itemp(i,j-1)+itemp(i+1,j-1) & |
|
992 -itemp(i-1,j+1)-2.*itemp(i,j+1)-itemp(i+1,j+1)) & |
|
993 /(8.0*30.0) |
|
994 ! |
|
995 y = (itemp(i-1,j-1)+2.*itemp(i-1,j)+itemp(i-1,j+1) & |
|
996 -itemp(i+1,j-1)-2.*itemp(i+1,j)-itemp(i+1,j+1)) & |
|
997 /(8.0*30.0) |
|
998 ! |
|
999 slope = 90. - atand(sqrt(x*x + y*y)) |
|
1000 if (x .eq. 0. .and. y .eq. 0.) then |
|
1001 aspect = 360. |
|
1002 else |
|
1003 aspect = atan2d(x,y) |
|
1004 end if |
|
1005 ! |
|
1006 c = sind(alt)*sind(slope) + cosd(alt)*cosd(slope)*cosd(az-aspect) |
|
1007 ! |
|
1008 if (c .le. 0.) then |
|
1009 idispvec(k)=0 |
|
1010 else |
|
1011 idispvec(k)=nint(c*100.) |
|
1012 end if |
|
1013 ! |
|
1014 ! flood with water if not painted |
|
1015 if (itemp(i,j) .le. iflood) idispvec(k) = 107 |
|
1016 ! |
|
1017 ! print *, i,j,k,idispvec(k) |
|
1018 k = k + 1 |
|
1019 end do |
|
1020 write (fmt,'("(",I0,"(i3)$)")') irowcols+1 |
|
1021 write (50,'(2A$)') ' "' |
|
1022 !write (50,108) (idispvec(l), l=0, irowcols) |
|
1023 write (50,fmt) (idispvec(l), l=0, irowcols) |
|
1024 write (50,'(2A)') '",' |
|
1025 !108 format (' "',<irowcols+1>(i3),'",') |
|
1026 k = 0 |
|
1027 end do |
|
1028 ! |
|
1029 write (50,112) |
|
1030 112 format ('};') |
|
1031 close (50) |
|
1032 ! |
|
1033 return |
|
1034 end subroutine XPMRelief |
|
1035 ! |
|
1036 !------------------------------------------------------------------------------ |
|
1037 ! |
|
1038 !********************************************************** |
|
1039 ! The following were added by Forrest Hoffman |
|
1040 ! Sat Apr 7 23:20:48 EDT 2007 |
|
1041 !********************************************************** |
|
1042 ! |
|
1043 !------------------------------------------------------------------------------ |
|
1044 ! |
|
1045 real function cosd(angle) |
|
1046 implicit none |
|
1047 real, intent(in) :: angle |
|
1048 |
|
1049 cosd = cos(angle * rad_per_deg) |
|
1050 |
|
1051 return |
|
1052 end function cosd |
|
1053 ! |
|
1054 !------------------------------------------------------------------------------ |
|
1055 ! |
|
1056 real function sind(angle) |
|
1057 implicit none |
|
1058 real, intent(in) :: angle |
|
1059 |
|
1060 sind = sin(angle * rad_per_deg) |
|
1061 |
|
1062 return |
|
1063 end function sind |
|
1064 ! |
|
1065 !------------------------------------------------------------------------------ |
|
1066 ! |
|
1067 real function atand(x) |
|
1068 implicit none |
|
1069 real, intent(in) :: x |
|
1070 |
|
1071 atand = atan(x) * deg_per_rad |
|
1072 |
|
1073 return |
|
1074 end function atand |
|
1075 ! |
|
1076 !------------------------------------------------------------------------------ |
|
1077 ! |
|
1078 real function atan2d(x, y) |
|
1079 implicit none |
|
1080 real, intent(in) :: x |
|
1081 real, intent(in) :: y |
|
1082 |
|
1083 atan2d = atan2(x, y) * deg_per_rad |
|
1084 |
|
1085 return |
|
1086 end function atan2d |
|
1087 ! |
|
1088 !------------------------------------------------------------------------------ |
|
1089 ! |
|
1090 |
|
1091 end module realizer_mod |