Context Navigation

module_generic.f90 @ 2341

Last change on this file since 2341 was 2340, checked in by lfita, 6 years ago

Adding:

`continguos_homogene_zones': Subroutine to look for contiguous zones by looking by continuous grid points
`get_xyconlimits': Subroutine for getting the limits of contiguous values from a given point in a 2D matrix
`Index1DArrayL': Function to provide the first index of a given value inside a 1D boolean array
`multi_Index1DArrayL': Subroutine to provide the indices of a given value inside a 1D boolean array

File size: 28.7 KB

Line
1	MODULE module_generic
2	! Module with generic functions
3
4	!!!!!!! Subroutines/Functions
5	! continguos_homogene_zones: Subroutine to look for contiguous zones by looking by continuous grid points
6	! freeunit: provides the number of a free unit in which open a file
7	! GetInNamelist: Subroutine to get a paramter from a namelistfile
8	! get_xyconlimits: Subroutine for getting the limits of contiguous values from a given point in a 2D matrix
9	! index_list_coordsI: Function to provide the index of a given coordinate within a list of integer coordinates
10	! Index1DArrayI: Function to provide the first index of a given value inside a 1D integer array
11	! Index1DArrayR: Function to provide the first index of a given value inside a 1D real array
12	! Index1DArrayR_K: Function to provide the first index of a given value inside a 1D real(r_k) array
13	! Index1DArrayL: Function to provide the first index of a given value inside a 1D boolean array
14	! Index2DArrayR: Function to provide the first index of a given value inside a 2D real array
15	! Index2DArrayR_K: Function to provide the first index of a given value inside a 2D real(r_k) array
16	! Nvalues_2DArrayI: Number of different values of a 2D integer array
17	! mat2DPosition: Function to provide the i, j indices of a given value inside a 2D matrix
18	! multi_Index1DArrayL: Subroutine to provide the indices of a given value inside a 1D boolean array
19	! numberTimes: Function to provide the number of times that a given set of characters happen within a string
20	! RangeI: Function to provide a range of d1 values from 'iniv' to 'endv', of integer values in a vector
21	! RangeR: Function to provide a range of d1 values from 'iniv' to 'endv', of real values in a vector
22	! RangeR_K: Function to provide a range of d1 from 'iniv' to 'endv', of real(r_k) values in a vector
23	! stoprun: Subroutine to stop running and print a message
24	! vectorI_S: Function to transform a vector of integers to a string of characters
25	! vectorR_S: Function to transform a vector of reals to a string of characters
26	! xzones_homogenization: Subroutine to homogenize 2D contiguous zones along x-axis. Zones might be
27	! contiguous, but with different number assigned !
28	! yzones_homogenization: Subroutine to homogenize 2D contiguous zones along y-axis. Zones might be
29	! contiguous, but with different number assigned !
30
31	USE module_definitions
32	USE module_basic
33
34	CONTAINS
35
36	SUBROUTINE Nvalues_2DArrayI(dx, dy, dxy, mat2DI, Nvals, vals)
37	! Subroutine to give the number of different values of a 2D integer array
38
39	IMPLICIT NONE
40
41	INTEGER, INTENT(in) :: dx, dy, dxy
42	INTEGER, DIMENSION(dx,dy), INTENT(in) :: mat2DI
43	INTEGER, INTENT(out) :: Nvals
44	INTEGER, DIMENSION(dxy), INTENT(out) :: vals
45
46	! Local
47	INTEGER :: i, j, ij
48
49	!!!!!!! Variables
50	! dx, dy: size of the 2D space
51	! mat2DI: 2D integer matrix
52	! Nvals: number of different values
53	! vals: vector with the different values
54
55	fname = 'Nvalues_2DArrayI'
56
57	vals = 0
58
59	Nvals = 1
60	vals(1) = mat2DI(1,1)
61	DO i=1,dx
62	DO j=1,dy
63	IF (Index1DArrayI(vals, Nvals, mat2DI(i,j)) == -1) THEN
64	Nvals = Nvals + 1
65	vals(Nvals) = mat2DI(i,j)
66	END IF
67	END DO
68	END DO
69
70	RETURN
71
72	END SUBROUTINE Nvalues_2DArrayI
73
74	INTEGER FUNCTION index_list_coordsI(Ncoords, coords, icoord)
75	! Function to provide the index of a given coordinate within a list of integer coordinates
76
77	IMPLICIT NONE
78
79	INTEGER, INTENT(in) :: Ncoords
80	INTEGER, DIMENSION(Ncoords,2), INTENT(in) :: coords
81	INTEGER, DIMENSION(2), INTENT(in) :: icoord
82
83	! Local
84	INTEGER, DIMENSION(Ncoords) :: dist
85	INTEGER :: i,mindist
86	INTEGER, DIMENSION(1) :: iloc
87
88	!!!!!!! Variables
89	! Ncoords: number of coordinates in the list
90	! coords: list of coordinates
91	! icoord: coordinate to find
92
93	fname = 'index_list_coordsI'
94
95	dist = (coords(:,1)-icoord(1))2+(coords(:,2)-icoord(2))2
96
97	IF (ANY(dist == 0)) THEN
98	iloc = MINLOC(dist)
99	index_list_coordsI = iloc(1)
100	ELSE
101	index_list_coordsI = -1
102	END IF
103
104	END FUNCTION index_list_coordsI
105
106	INTEGER FUNCTION Index1DArrayI(array1D, d1, val)
107	! Function to provide the first index of a given value inside a 1D integer array
108
109	IMPLICIT NONE
110
111	INTEGER, INTENT(in) :: d1
112	INTEGER, INTENT(in) :: val
113	INTEGER, DIMENSION(d1), INTENT(in) :: array1D
114
115	! Local
116	INTEGER :: i
117
118	fname = 'Index1DArrayI'
119
120	Index1DArrayI = -1
121
122	DO i=1,d1
123	IF (array1d(i) == val) THEN
124	Index1DArrayI = i
125	EXIT
126	END IF
127	END DO
128
129	END FUNCTION Index1DArrayI
130
131	INTEGER FUNCTION Index1DArrayR(array1D, d1, val)
132	! Function to provide the first index of a given value inside a 1D real array
133
134	IMPLICIT NONE
135
136	INTEGER, INTENT(in) :: d1
137	REAL, INTENT(in) :: val
138	REAL, DIMENSION(d1), INTENT(in) :: array1D
139
140	! Local
141	INTEGER :: i
142
143	fname = 'Index1DArrayR'
144
145	Index1DArrayR = -1
146
147	DO i=1,d1
148	IF (array1d(i) == val) THEN
149	Index1DArrayR = i
150	EXIT
151	END IF
152	END DO
153
154	END FUNCTION Index1DArrayR
155
156	INTEGER FUNCTION Index1DArrayR_K(array1D, d1, val)
157	! Function to provide the first index of a given value inside a 1D real(r_k) array
158
159	IMPLICIT NONE
160
161	INTEGER, INTENT(in) :: d1
162	REAL(r_k), INTENT(in) :: val
163	REAL(r_k), DIMENSION(d1), INTENT(in) :: array1D
164
165	! Local
166	INTEGER :: i
167
168	fname = 'Index1DArrayR_K'
169
170	Index1DArrayR_K = -1
171
172	DO i=1,d1
173	IF (array1d(i) == val) THEN
174	Index1DArrayR_K = i
175	EXIT
176	END IF
177	END DO
178
179	END FUNCTION Index1DArrayR_K
180
181	INTEGER FUNCTION Index1DArrayL(array1D, d1, val)
182	! Function to provide the first index of a given value inside a 1D boolean array
183
184	IMPLICIT NONE
185
186	INTEGER, INTENT(in) :: d1
187	LOGICAL, INTENT(in) :: val
188	LOGICAL, DIMENSION(d1), INTENT(in) :: array1D
189
190	! Local
191	INTEGER :: i
192	CHARACTER(LEN=50) :: fname
193
194	fname = 'Index1DArrayL'
195
196	Index1DArrayL = -1
197
198	DO i=1,d1
199	IF (array1d(i) .EQV. val) THEN
200	Index1DArrayL = i
201	EXIT
202	END IF
203	END DO
204
205	END FUNCTION Index1DArrayL
206
207	SUBROUTINE multi_Index1DArrayL(array1D, d1, val, Ntimes, pos)
208	! Subroutine to provide the indices of a given value inside a 1D boolean array
209
210	IMPLICIT NONE
211
212	INTEGER, INTENT(in) :: d1
213	LOGICAL, INTENT(in) :: val
214	LOGICAL, DIMENSION(d1), INTENT(in) :: array1D
215	INTEGER, INTENT(out) :: Ntimes
216	INTEGER, DIMENSION(d1), INTENT(out) :: pos
217
218	! Local
219	INTEGER :: i
220	CHARACTER(LEN=50) :: fname
221
222	!!!!!!! Variables
223	! array1D: 1D Array of values
224	! d1: length of array
225	! val: Value to look for
226	! Ntimes: Number of times val is found within array1D
227	! pos: positions of the values of val
228
229	fname = 'multi_Index1DArrayL'
230
231	Ntimes = 0
232	pos = -1
233
234	DO i=1,d1
235	IF (array1d(i) .EQV. val) THEN
236	Ntimes = Ntimes + 1
237	pos(Ntimes) = i
238	END IF
239	END DO
240
241	END SUBROUTINE multi_Index1DArrayL
242
243	FUNCTION Index2DArrayR(array2D, d1, d2, val)
244	! Function to provide the first index of a given value inside a 2D real array
245
246	IMPLICIT NONE
247
248	INTEGER, INTENT(in) :: d1, d2
249	REAL, INTENT(in) :: val
250	REAL, DIMENSION(d1,d2), INTENT(in) :: array2D
251	INTEGER, DIMENSION(2) :: Index2DArrayR
252
253	! Local
254	INTEGER :: i, j
255
256	fname = 'Index2DArrayR'
257
258	Index2DArrayR = -1
259
260	DO i=1,d1
261	DO j=1,d2
262	IF (array2d(i,j) == val) THEN
263	Index2DArrayR(1) = i
264	Index2DArrayR(2) = j
265	EXIT
266	END IF
267	END DO
268	END DO
269
270	END FUNCTION Index2DArrayR
271
272	FUNCTION Index2DArrayR_K(array2D, d1, d2, val)
273	! Function to provide the first index of a given value inside a 2D real array
274
275	IMPLICIT NONE
276
277	INTEGER, INTENT(in) :: d1, d2
278	REAL(r_k), INTENT(in) :: val
279	REAL(r_k), DIMENSION(d1,d2), INTENT(in) :: array2D
280	INTEGER, DIMENSION(2) :: Index2DArrayR_K
281
282	! Local
283	INTEGER :: i, j
284
285	fname = 'Index2DArrayR_K'
286
287	Index2DArrayR_K = -1
288
289	DO i=1,d1
290	DO j=1,d2
291	IF (array2d(i,j) == val) THEN
292	Index2DArrayR_K(1) = i
293	Index2DArrayR_K(2) = j
294	EXIT
295	END IF
296	END DO
297	END DO
298
299	END FUNCTION Index2DArrayR_K
300
301	INTEGER FUNCTION numberTimes(String, charv)
302	! Function to provide the number of times that a given set of characters happen within a string
303
304	IMPLICIT NONE
305
306	CHARACTER(LEN=*), INTENT(IN) :: String, charv
307
308	! Local
309	INTEGER :: i, Lstring, Lcharv
310
311	numberTimes = 0
312
313	Lstring = LEN_TRIM(String)
314	Lcharv = LEN_TRIM(charv)
315
316	DO i=1,Lstring - Lcharv
317	IF (String(i:i+Lcharv-1) == TRIM(charv)) numberTimes = numberTimes + 1
318	END DO
319
320	END FUNCTION numberTimes
321
322
323	FUNCTION RangeI(d1, iniv, endv)
324	! Function to provide a range of d1 values from 'iniv' to 'endv', of integer values in a vector
325
326	IMPLICIT NONE
327
328	INTEGER, INTENT(in) :: d1, iniv, endv
329	INTEGER, DIMENSION(d1) :: RangeI
330
331	! Local
332	INTEGER :: i, intv
333
334	fname = 'RangeI'
335
336	intv = (endv - iniv) / (d1*1 - 1)
337
338	RangeI(1) = iniv
339	DO i=2,d1
340	RangeI(i) = RangeI(i-1) + intv
341	END DO
342
343	END FUNCTION RangeI
344
345	FUNCTION RangeR(d1, iniv, endv)
346	! Function to provide a range of d1 from 'iniv' to 'endv', of real values in a vector
347
348	IMPLICIT NONE
349
350	INTEGER, INTENT(in) :: d1
351	REAL, INTENT(in) :: iniv, endv
352	REAL, DIMENSION(d1) :: RangeR
353
354	! Local
355	INTEGER :: i
356	REAL :: intv
357
358	fname = 'RangeR'
359
360	intv = (endv - iniv) / (d1*1. - 1.)
361
362	RangeR(1) = iniv
363	DO i=2,d1
364	RangeR(i) = RangeR(i-1) + intv
365	END DO
366
367	END FUNCTION RangeR
368
369	FUNCTION RangeR_K(d1, iniv, endv)
370	! Function to provide a range of d1 from 'iniv' to 'endv', of real(r_k) values in a vector
371
372	IMPLICIT NONE
373
374	INTEGER, INTENT(in) :: d1
375	REAL(r_k), INTENT(in) :: iniv, endv
376	REAL(r_k), DIMENSION(d1) :: RangeR_K
377
378	! Local
379	INTEGER :: i
380	REAL(r_k) :: intv
381
382	fname = 'RangeR_K'
383
384	intv = (endv - iniv) / (d1*oneRK-oneRK)
385
386	RangeR_K(1) = iniv
387	DO i=2,d1
388	RangeR_K(i) = RangeR_K(i-1) + intv
389	END DO
390
391	END FUNCTION RangeR_K
392
393	INTEGER FUNCTION freeunit()
394	! provides the number of a free unit in which open a file
395
396	IMPLICIT NONE
397
398	LOGICAL :: is_used
399
400	is_used = .true.
401	DO freeunit=10,100
402	INQUIRE(unit=freeunit, opened=is_used)
403	IF (.not. is_used) EXIT
404	END DO
405
406	RETURN
407
408	END FUNCTION freeunit
409
410	SUBROUTINE GetInNamelist(namelistfile, param, kindparam, Ival, Rval, Lval, Sval)
411	! Subroutine to get a paramter from a namelistfile
412
413	IMPLICIT NONE
414
415	CHARACTER(LEN=*), INTENT(IN) :: namelistfile, param
416	CHARACTER(LEN=1), INTENT(IN) :: kindparam
417	INTEGER, OPTIONAL, INTENT(OUT) :: Ival
418	REAL, OPTIONAL, INTENT(OUT) :: Rval
419	LOGICAL, OPTIONAL, INTENT(OUT) :: Lval
420	CHARACTER(LEN=200), OPTIONAL, INTENT(OUT) :: Sval
421
422	! Local
423	INTEGER :: i, funit, ios
424	INTEGER :: Lparam, posparam
425	LOGICAL :: is_used
426	CHARACTER(LEN=1000) :: line, message
427	CHARACTER(LEN=200), DIMENSION(2) :: lvals
428	CHARACTER(LEN=200) :: pval
429
430	!!!!!!! Variables
431	! namelistfile: name of the namelist file
432	! param: parameter to get
433	! paramkind: kind of the parameter (I: Integer, L: boolean, R: Real, S: String)
434
435	fname = 'GetInNamelist'
436
437	! Reading dimensions file and defining dimensions
438	is_used = .true.
439	DO funit=10,100
440	INQUIRE(unit=funit, opened=is_used)
441	IF (.not. is_used) EXIT
442	END DO
443
444	OPEN(funit, FILE=TRIM(namelistfile), STATUS='old', FORM='formatted', IOSTAT=ios)
445	IF ( ios /= 0 ) CALL stoprun(message, fname)
446
447	Lparam = LEN_TRIM(param)
448
449	DO i=1,10000
450	READ(funit,"(A200)",END=100)line
451	posparam = INDEX(TRIM(line), TRIM(param))
452	IF (posparam /= 0) EXIT
453
454	END DO
455	100 CONTINUE
456
457	IF (posparam == 0) THEN
458	message = "namelist '" // TRIM(namelistfile) // "' does not have parameter '" // TRIM(param) // &
459	"' !!"
460	CALL stoprun(message, fname)
461	END IF
462
463	CLOSE(UNIT=funit)
464
465	CALL split(line, '=', 2, lvals)
466	IF (kindparam /= 'S') THEN
467	CALL RemoveNonNum(lvals(2), pval)
468	END IF
469
470	! L. Fita, LMD. October 2015
471	! Up to now, only getting scalar values
472	kparam: SELECT CASE (kindparam)
473	CASE ('I')
474	Ival = StoI(pval)
475	! PRINT *,TRIM(param),'= ', Ival
476	CASE ('L')
477	Lval = StoL(pval)
478	! PRINT *,TRIM(param),'= ', Lval
479	CASE ('R')
480	Rval = StoR(pval)
481	! PRINT *,TRIM(param),'= ', Rval
482	CASE ('S')
483	Sval = lvals(2)
484
485	CASE DEFAULT
486	message = "type of parameter '" // kindparam // "' not ready !!"
487	CALL stoprun(message, fname)
488
489	END SELECT kparam
490
491	END SUBROUTINE GetInNamelist
492
493	SUBROUTINE stoprun(msg, fname)
494	! Subroutine to stop running and print a message
495
496	IMPLICIT NONE
497
498	CHARACTER(LEN=*), INTENT(IN) :: fname
499	CHARACTER(LEN=*), INTENT(IN) :: msg
500
501	! local
502	CHARACTER(LEN=50) :: errmsg, warnmsg
503
504	errmsg = 'ERROR -- error -- ERROR -- error'
505
506	PRINT *, TRIM(errmsg)
507	PRINT *, ' ' // TRIM(fname) // ': ' // TRIM(msg)
508	STOP
509
510	END SUBROUTINE stoprun
511
512	SUBROUTINE xzones_homogenization(dx, dy, inzones, outzones)
513	! Subroutine to homogenize 2D contiguous zones along x-axis, zones might be contiguous, but with
514	! different number assigned !
515	! Here we have a 2D matrix of integers, with contiguous integer filled zones, zero outside any zone
516	! It might be, that within the same zone, might be lines which do not share the same integer
517	! 0 0 0 0 0 0 0 0 0 0
518	! 0 1 1 0 0 0 1 1 0 0
519	! 0 2 0 0 1 == > 0 1 0 0 2
520	! 0 1 1 0 0 0 1 1 0 0
521
522	IMPLICIT NONE
523
524	INTEGER, INTENT(in) :: dx, dy
525	INTEGER, DIMENSION(dx,dy), INTENT(in) :: inzones
526	INTEGER, DIMENSION(dx,dy), INTENT(out) :: outzones
527
528	! Local
529	INTEGER :: i,j,k
530	INTEGER :: Nmaxzones, TOTzones
531	LOGICAL :: assigned
532	INTEGER, DIMENSION(dx) :: prevline
533	INTEGER, DIMENSION(dy) :: Nxzones
534	INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: zones
535
536	!!!!!!! Variables
537	! dx, dy: Shape of the 2D space
538	! inzones: zones to homogenize
539	! outzones: zones homogenized
540
541	fname = 'xzones_homogenization'
542
543	! Maximum possible number of zones
544	Nmaxzones = INT((dx/2)*(dy/2))
545
546	! Matrix with [i,j,Nzone,izone/ezone]
547	IF (ALLOCATED(zones)) DEALLOCATE(zones)
548	ALLOCATE(zones(dy,Nmaxzones,3))
549
550	zones = 0
551	Nxzones = 0
552	! Getting beginning/end of y-bands
553	DO j=1, dy
554	k = 0
555	i = 1
556	IF (inzones(i,j) /= 0) THEN
557	k = k + 1
558	zones(j,k,1) = i
559	zones(j,k,3) = k
560	END IF
561	DO i=2, dx
562	IF ( (inzones(i,j) /= 0) .AND. (inzones(i-1,j) == 0)) THEN
563	k = k+1
564	zones(j,k,1) = i
565	zones(j,k,3) = k
566	ELSE IF ( (inzones(i-1,j) /= 0) .AND. (inzones(i,j) == 0)) THEN
567	zones(j,k,2) = i-1
568	zones(j,k,3) = k
569	END IF
570	END DO
571	IF (k > 0) THEN
572	IF (zones(j,k,2) == 0) zones(j,k,2) = dx
573	END IF
574	Nxzones(j) = k
575	END DO
576
577	! Homogenizing contigous zones
578	outzones = 0
579	TOTzones = 0
580	j = 1
581	DO k = 1, Nxzones(j)
582	TOTzones = TOTzones + 1
583	DO i=zones(j,k,1), zones(j,k,2)
584	outzones(i,j) = TOTzones
585	END DO
586	END DO
587
588	DO j=2, dy
589	prevline = outzones(:,j-1)
590	DO k = 1, Nxzones(j)
591	assigned = .FALSE.
592	DO i=zones(j,k,1), zones(j,k,2)
593	IF (prevline(i) /= 0) THEN
594	outzones(zones(j,k,1):zones(j,k,2),j) = prevline(i)
595	assigned = .TRUE.
596	EXIT
597	END IF
598	END DO
599	IF (.NOT.assigned) THEN
600	TOTzones = TOTzones + 1
601	DO i=zones(j,k,1), zones(j,k,2)
602	outzones(i,j) = TOTzones
603	END DO
604	END IF
605	END DO
606	END DO
607
608	IF (ALLOCATED(zones)) DEALLOCATE(zones)
609
610	END SUBROUTINE xzones_homogenization
611
612	SUBROUTINE yzones_homogenization(dx, dy, inzones, outzones)
613	! Subroutine to homogenize 2D contiguous zones along y-axis, zones might be contiguous, but with
614	! different number assigned !
615	! Here we have a 2D matrix of integers, with contiguous integer filled zones, zero outside any zone
616	! It might be, that within the same zone, might be lines which do not share the same integer
617	! 0 0 0 0 0 0 0 0 0 0
618	! 0 1 1 0 0 0 1 1 0 0
619	! 0 2 0 0 1 == > 0 1 0 0 2
620	! 0 1 1 0 0 0 1 1 0 0
621
622	IMPLICIT NONE
623
624	INTEGER, INTENT(in) :: dx, dy
625	INTEGER, DIMENSION(dx,dy), INTENT(in) :: inzones
626	INTEGER, DIMENSION(dx,dy), INTENT(out) :: outzones
627
628	! Local
629	INTEGER :: i,j,k
630	INTEGER :: Nmaxzones, TOTzones
631	LOGICAL :: assigned
632	INTEGER, DIMENSION(dy) :: prevline
633	INTEGER, DIMENSION(dx) :: Nyzones
634	INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: zones
635
636	!!!!!!! Variables
637	! dx, dy: Shape of the 2D space
638	! inzones: zones to homogenize
639	! outzones: zones homogenized
640
641	fname = 'yzones_homogenization'
642
643	! Maximum possible number of zones
644	Nmaxzones = INT((dx/2)*(dy/2))
645
646	! Matrix with [i,j,Nzone,izone/ezone]
647	IF (ALLOCATED(zones)) DEALLOCATE(zones)
648	ALLOCATE(zones(dx,Nmaxzones,3))
649
650	zones = 0
651	Nyzones = 0
652	! Getting beginning/end of y-bands
653	DO i=1, dx
654	k = 0
655	j = 1
656	IF (inzones(i,j) /= 0) THEN
657	k = k + 1
658	zones(i,k,1) = j
659	zones(i,k,3) = k
660	END IF
661	DO j=2, dy
662	IF ( (inzones(i,j) /= 0) .AND. (inzones(i,j-1) == 0)) THEN
663	k = k+1
664	zones(i,k,1) = j
665	zones(i,k,3) = k
666	ELSE IF ( (inzones(i,j-1) /= 0) .AND. (inzones(i,j) == 0)) THEN
667	zones(i,k,2) = j-1
668	zones(i,k,3) = k
669	END IF
670	END DO
671	IF (k > 0) THEN
672	IF (zones(i,k,2) == 0) zones(i,k,2) = dy
673	END IF
674	Nyzones(i) = k
675	END DO
676
677	! Homogenizing contigous zones
678	outzones = 0
679	TOTzones = 0
680	i = 1
681	DO k = 1, Nyzones(i)
682	TOTzones = TOTzones + 1
683	DO j=zones(i,k,1), zones(i,k,2)
684	outzones(i,j) = TOTzones
685	END DO
686	END DO
687
688	DO i=2, dx
689	prevline = outzones(i-1,:)
690	DO k = 1, Nyzones(i)
691	assigned = .FALSE.
692	DO j=zones(i,k,1), zones(i,k,2)
693	IF (prevline(j) /= 0) THEN
694	outzones(i,zones(i,k,1):zones(i,k,2)) = prevline(j)
695	assigned = .TRUE.
696	EXIT
697	END IF
698	END DO
699	IF (.NOT.assigned) THEN
700	TOTzones = TOTzones + 1
701	DO j=zones(i,k,1), zones(i,k,2)
702	outzones(i,j) = TOTzones
703	END DO
704	END IF
705	END DO
706	END DO
707
708	IF (ALLOCATED(zones)) DEALLOCATE(zones)
709
710	END SUBROUTINE yzones_homogenization
711
712	CHARACTER(len=1000) FUNCTION vectorI_S(d1, vector)
713	! Function to transform a vector of integers to a string of characters
714
715	IMPLICIT NONE
716
717	INTEGER, INTENT(in) :: d1
718	INTEGER, DIMENSION(d1), INTENT(in) :: vector
719
720	! Local
721	INTEGER :: iv
722	CHARACTER(len=50) :: IS
723
724	!!!!!!! Variables
725	! d1: length of the vector
726	! vector: values to transform
727
728	fname = 'vectorI_S'
729
730	vectorI_S = ''
731	DO iv=1, d1
732	WRITE(IS, '(I50)')vector(iv)
733	IF (iv == 1) THEN
734	vectorI_S = TRIM(IS)
735	ELSE
736	vectorI_S = TRIM(vectorI_S) // ', ' // TRIM(IS)
737	END IF
738	END DO
739
740	END FUNCTION vectorI_S
741
742	CHARACTER(len=1000) FUNCTION vectorR_S(d1, vector)
743	! Function to transform a vector of reals to a string of characters
744
745	IMPLICIT NONE
746
747	INTEGER, INTENT(in) :: d1
748	REAL, DIMENSION(d1), INTENT(in) :: vector
749
750	! Local
751	INTEGER :: iv
752	CHARACTER(len=50) :: RS
753
754	!!!!!!! Variables
755	! d1: length of the vector
756	! vector: values to transform
757
758	fname = 'vectorR_S'
759
760	vectorR_S = ''
761	DO iv=1, d1
762	WRITE(RS, '(F50.25)')vector(iv)
763	IF (iv == 1) THEN
764	vectorR_S = TRIM(RS)
765	ELSE
766	vectorR_S = TRIM(vectorR_S) // ', ' // TRIM(RS)
767	END IF
768	END DO
769
770	END FUNCTION vectorR_S
771
772	SUBROUTINE continguos_homogene_zones(dx, dy, matvals, Nzones, contzones)
773	! Subroutine to look for contiguous zones by looking by continuous grid points
774
775	IMPLICIT NONE
776
777	INTEGER, INTENT(in) :: dx, dy
778	INTEGER, DIMENSION(dx,dy), INTENT(in) :: matvals
779	INTEGER, INTENT(out) :: Nzones
780	INTEGER, DIMENSION(dx,dy), INTENT(out) :: contzones
781	! Local
782	INTEGER :: i,j, k
783	INTEGER :: ii, ei, ij, ej
784	INTEGER :: Ncont, Nassigned, Ncontmin
785	LOGICAL, DIMENSION(dx,dy) :: assigned, notdone
786	INTEGER, DIMENSION(:), ALLOCATABLE :: pzones, allzones
787	LOGICAL, DIMENSION(:), ALLOCATABLE :: passigns
788
789	!!!!!!! Variables
790	! dx, dy: shape of the matrix
791	! matvals: matrix with the values
792	! contzones: homogeneous zones found
793
794	fname = 'continguos_homogene_zones'
795
796	! Vector to keep track of all zone values
797	IF (ALLOCATED(allzones)) DEALLOCATE(allzones)
798	ALLOCATE(allzones(dx*dy/4))
799	allzones = 0
800
801	assigned = .FALSE.
802	notdone = .TRUE.
803	contzones = -1
804	Nzones = 0
805	DO i=1, dx
806	DO j=1, dy
807	! First
808	IF (matvals(i,j) /= 0 .AND. notdone(i,j)) THEN
809	CALL get_xyconlimits(dx, dy, matvals, 0, i, j, ii, ei, ij, ej)
810
811	! Has any point of the rays already been assigned?
812	! Along x
813	IF (ALLOCATED(pzones)) DEALLOCATE(pzones)
814	ALLOCATE(pzones(ei-ii+1))
815	IF (ALLOCATED(passigns)) DEALLOCATE(passigns)
816	ALLOCATE(passigns(ei-ii+1))
817	passigns = assigned(ii:ei,j)
818	CALL multi_Index1DArrayL(passigns, ei-ii+1, .TRUE., Nassigned, pzones)
819	IF (Nassigned /= 0) THEN
820	Ncontmin = 10000000
821	DO k=1, Nassigned
822	IF (contzones(ii+pzones(k)-1,j) < Ncontmin) Ncontmin = contzones(ii+pzones(k)-1,j)
823	END DO
824	! If there is more than one assigned value change all values to the minimum one
825	DO k=1, Nassigned
826	IF (contzones(ii+pzones(k)-1,j) /= Ncontmin) THEN
827	WHERE (contzones == contzones(ii+pzones(k)-1,j) .AND. contzones /= Ncontmin)
828	contzones = Ncontmin
829	END WHERE
830	allzones(contzones(ii+pzones(k)-1,j)) = 0
831	END IF
832	END DO
833	Ncont = Ncontmin
834	END IF
835	! Along y
836	IF (ALLOCATED(pzones)) DEALLOCATE(pzones)
837	ALLOCATE(pzones(ej-ij+1))
838	IF (ALLOCATED(passigns)) DEALLOCATE(passigns)
839	ALLOCATE(passigns(ej-ij+1))
840	passigns = assigned(i,ij:ej)
841	CALL multi_Index1DArrayL(passigns, ej-ij+1, .TRUE., Nassigned, pzones)
842	IF (Nassigned /= 0) THEN
843	Ncontmin = 10000000
844	DO k=1, Nassigned
845	IF (contzones(i,ij+pzones(k)-1) < Ncontmin) Ncontmin = contzones(i,ij+pzones(k)-1)
846	END DO
847	! If there is more than one assigned value change all values to the minimum one
848	DO k=1, Nassigned
849	IF (contzones(i,ij+pzones(k)-1) /= Ncontmin) THEN
850	WHERE (contzones == contzones(i,ij+pzones(k)-1) .AND. contzones /= Ncontmin)
851	contzones = Ncontmin
852	END WHERE
853	allzones(contzones(i,ij+pzones(k)-1)) = 0
854	END IF
855	END DO
856	Ncont = Ncontmin
857	END IF
858
859	IF (.NOT.assigned(i,j)) THEN
860	Nzones = Nzones + 1
861	Ncont = Nzones
862	allzones(Nzones) = 1
863	ELSE
864	Ncont = contzones(i,j)
865	END IF
866	contzones(i,j) = Ncont
867	contzones(ii:ei,j) = Ncont
868	contzones(i,ij:ej) = Ncont
869	notdone(i,j) = .FALSE.
870	assigned(i,j) = .TRUE.
871	assigned(ii:ei,j) = .TRUE.
872	assigned(i,ij:ej) = .TRUE.
873	END IF
874	END DO
875	END DO
876
877	! Using allzones to provide continuous assigned values
878	Nzones = 0
879	DO k=1, dx*dy/4
880	IF (allzones(k) /= 0) THEN
881	Nzones = Nzones + 1
882	IF (allzones(k) /= Nzones) THEN
883	WHERE(contzones == allzones(k))
884	contzones = Nzones
885	END WHERE
886	END IF
887	END IF
888	END DO
889
890	RETURN
891
892	END SUBROUTINE continguos_homogene_zones
893
894	SUBROUTINE get_xyconlimits(d1, d2, matv, NOval, i, j, ix, ex, iy, ey)
895	! Subroutine for getting the limits of contiguous values from a given point in a 2D matrix
896
897	IMPLICIT NONE
898
899	INTEGER, INTENT(in) :: d1, d2, i, j, NOval
900	INTEGER, DIMENSION(d1,d2), INTENT(in) :: matv
901	INTEGER, INTENT(out) :: ix, ex, iy, ey
902
903	! Local
904	INTEGER :: i1, j1
905	LOGICAL :: found
906
907	!!!!!!! Variables
908	! d1, d2: Shape of input 2D values
909	! i, j: grid point from which get the contiguous values
910	! NOval: value for grid points without data
911	! matv: 2D matrx values
912	! ix, ex, iy, ey: limits from a grid point
913
914	fname = 'get_xyconlimits'
915
916	ix = -1
917	ex = -1
918	iy = -1
919	ey = -1
920
921	! Before i
922	IF (i > 1) THEN
923	ix = i
924	found = .FALSE.
925	DO i1=i-1,1,-1
926	IF (matv(i1,j) == NOval) THEN
927	ix = i1 + 1
928	found = .TRUE.
929	EXIT
930	END IF
931	END DO
932	IF (.NOT.found) ix=1
933	ELSE
934	ix = 1
935	END IF
936	! After i
937	IF (i < d1) THEN
938	ex = i
939	found = .FALSE.
940	DO i1=i+1,d1
941	IF (matv(i1,j) == NOval) THEN
942	ex = i1 - 1
943	found = .TRUE.
944	EXIT
945	END IF
946	END DO
947	IF (.NOT.found) ex=d1
948	ELSE
949	ex = d1
950	END IF
951
952	! Before j
953	IF (j > 1) THEN
954	iy = j
955	found = .FALSE.
956	DO j1=j-1,1,-1
957	IF (matv(i,j1) == NOval) THEN
958	iy = j1 + 1
959	found = .TRUE.
960	EXIT
961	END IF
962	END DO
963	IF (.NOT.found) iy=1
964	ELSE
965	iy = 1
966	END IF
967	! After j
968	IF (j < d2) THEN
969	ey = j
970	found = .FALSE.
971	DO j1=j+1,d2
972	IF (matv(i,j1) == NOval) THEN
973	ey = j1 - 1
974	found = .TRUE.
975	EXIT
976	END IF
977	END DO
978	IF (.NOT.found) ey=d2
979	ELSE
980	ey = d2
981	END IF
982
983	RETURN
984
985	END SUBROUTINE get_xyconlimits
986
987	END MODULE module_generic

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format