001 /*
002 * $Id: ProfileAlongTrack.java,v 1.12 2012/02/19 17:35:42 davep Exp $
003 *
004 * This file is part of McIDAS-V
005 *
006 * Copyright 2007-2012
007 * Space Science and Engineering Center (SSEC)
008 * University of Wisconsin - Madison
009 * 1225 W. Dayton Street, Madison, WI 53706, USA
010 * https://www.ssec.wisc.edu/mcidas
011 *
012 * All Rights Reserved
013 *
014 * McIDAS-V is built on Unidata's IDV and SSEC's VisAD libraries, and
015 * some McIDAS-V source code is based on IDV and VisAD source code.
016 *
017 * McIDAS-V is free software; you can redistribute it and/or modify
018 * it under the terms of the GNU Lesser Public License as published by
019 * the Free Software Foundation; either version 3 of the License, or
020 * (at your option) any later version.
021 *
022 * McIDAS-V is distributed in the hope that it will be useful,
023 * but WITHOUT ANY WARRANTY; without even the implied warranty of
024 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
025 * GNU Lesser Public License for more details.
026 *
027 * You should have received a copy of the GNU Lesser Public License
028 * along with this program. If not, see http://www.gnu.org/licenses.
029 */
030
031 package edu.wisc.ssec.mcidasv.data.hydra;
032
033 import visad.Data;
034 import visad.FlatField;
035 import visad.Set;
036 import visad.Gridded1DSet;
037 import visad.CoordinateSystem;
038 import visad.RealType;
039 import visad.RealTupleType;
040 import visad.SetType;
041 import visad.Linear1DSet;
042 import visad.Linear2DSet;
043 import visad.Unit;
044 import visad.FunctionType;
045 import visad.VisADException;
046 import visad.QuickSort;
047 import java.rmi.RemoteException;
048
049 import java.util.Hashtable;
050 import java.util.HashMap;
051 import java.util.StringTokenizer;
052
053 import java.io.BufferedReader;
054 import java.io.FileInputStream;
055 import java.io.IOException;
056 import java.io.InputStream;
057 import java.io.InputStreamReader;
058
059
060 public abstract class ProfileAlongTrack extends MultiDimensionAdapter {
061
062 private FunctionType mathtype;
063
064 int TrackLen;
065 int VertLen;
066
067 private float[] vertLocs = null;
068 private float[] trackTimes = null;
069 private float[] trackLongitude = null;
070 private float[] trackLatitude = null;
071
072 public static String longitude_name = "Longitude";
073 public static String latitude_name = "Latitude";
074 public static String trackDim_name = "TrackDim";
075 public static String vertDim_name = "VertDim";
076 public static String array_name = "array_name";
077 public static String profileTime_name = "ProfileTime";
078 public static String profileTime_unit = "ProfileTime_unit";
079 public static String altitude_unit = "altitude_unit";
080 public static String sfcElev_name = "SurfaceElev";
081 public static String range_name = "range_name";
082 public static String scale_name = "scale_name";
083 public static String offset_name = "offset_name";
084 public static String fill_value_name = "fill_value_name";
085 public static String valid_range = "valid_range";
086 public static String ancillary_file_name = "ancillary_file";
087 static String product_name = "product_name";
088
089 String[] rangeName_s = null;
090 Class[] arrayType_s = null;
091 Unit[] rangeUnit_s = new Unit[] {null};
092
093 RealType track = RealType.getRealType(trackDim_name);
094 RealType vert = RealType.getRealType(vertDim_name);
095 RealType[] domainRealTypes = new RealType[2];
096
097 RealType vertLocType;
098 RealType trackTimeType;
099
100 int track_idx = -1;
101 int vert_idx = -1;
102 int range_rank = -1;
103
104 int track_tup_idx;
105 int vert_tup_idx;
106
107 boolean isVertDimAlt = true;
108
109 CoordinateSystem cs = null;
110
111 public static HashMap getEmptySubset() {
112 HashMap<String, double[]> subset = new HashMap<String, double[]>();
113 subset.put(trackDim_name, new double[3]);
114 subset.put(vertDim_name, new double[3]);
115 return subset;
116 }
117
118 public static HashMap getEmptyMetadataTable() {
119 HashMap<String, String> metadata = new HashMap<String, String>();
120 metadata.put(array_name, null);
121 metadata.put(trackDim_name, null);
122 metadata.put(vertDim_name, null);
123 metadata.put(longitude_name, null);
124 metadata.put(latitude_name, null);
125 metadata.put(profileTime_name, null);
126 metadata.put(profileTime_unit, null);
127 metadata.put(altitude_unit, null);
128 metadata.put(sfcElev_name, null);
129 metadata.put(scale_name, null);
130 metadata.put(offset_name, null);
131 metadata.put(fill_value_name, null);
132 /*
133 metadata.put(range_name, null);
134 metadata.put(range_unit, null);
135 metadata.put(valid_range, null);
136 */
137 return metadata;
138 }
139
140 public ProfileAlongTrack() {
141 }
142
143 public ProfileAlongTrack(MultiDimensionReader reader, HashMap metadata) {
144 this(reader, metadata, true);
145 }
146
147 public ProfileAlongTrack(MultiDimensionReader reader, HashMap metadata, boolean isVertDimAlt) {
148 super(reader, metadata);
149 this.isVertDimAlt = isVertDimAlt;
150 this.init();
151 }
152
153
154 private void init() {
155 for (int k=0; k<array_rank;k++) {
156 if ( ((String)metadata.get(trackDim_name)).equals(array_dim_names[k]) ) {
157 track_idx = k;
158 }
159 if ( ((String)metadata.get(vertDim_name)).equals(array_dim_names[k]) ) {
160 vert_idx = k;
161 }
162 }
163
164 int[] lengths = new int[2];
165
166 if (track_idx < vert_idx) {
167 domainRealTypes[0] = vert;
168 domainRealTypes[1] = track;
169 track_tup_idx = 1;
170 vert_tup_idx = 0;
171 lengths[0] = array_dim_lengths[vert_idx];
172 lengths[1] = array_dim_lengths[track_idx];
173 }
174 else {
175 domainRealTypes[0] = track;
176 domainRealTypes[1] = vert;
177 track_tup_idx = 0;
178 vert_tup_idx = 1;
179 lengths[0] = array_dim_lengths[track_idx];
180 lengths[1] = array_dim_lengths[vert_idx];
181 }
182
183 TrackLen = array_dim_lengths[track_idx];
184 VertLen = array_dim_lengths[vert_idx];
185
186 String rangeName = null;
187 if (metadata.get("range_name") != null) {
188 rangeName = (String)metadata.get("range_name");
189 }
190 else {
191 rangeName = (String)metadata.get("array_name");
192 }
193 rangeType = RealType.getRealType(rangeName, rangeUnit_s[0]);
194
195 try {
196 rangeProcessor = RangeProcessor.createRangeProcessor(reader, metadata);
197 }
198 catch (Exception e) {
199 System.out.println("RangeProcessor failed to create");
200 e.printStackTrace();
201 }
202
203 try {
204 if (isVertDimAlt) {
205 vertLocs = getVertBinAltitude();
206 }
207 vertLocType = makeVertLocType();
208 trackTimes = getTrackTimes();
209 trackTimeType = makeTrackTimeType();
210 trackLongitude = getTrackLongitude();
211 trackLatitude = getTrackLatitude();
212 }
213 catch (Exception e) {
214 System.out.println(e);
215 }
216
217 }
218
219 public int getTrackLength() {
220 return TrackLen;
221 }
222
223 public int getVertLength() {
224 return VertLen;
225 }
226
227 public int getVertIdx() {
228 return vert_idx;
229 }
230
231 public int getTrackIdx() {
232 return track_idx;
233 }
234
235 public int getVertTupIdx() {
236 return vert_tup_idx;
237 }
238
239 public int getTrackTupIdx() {
240 return track_tup_idx;
241 }
242
243 public Set makeDomain(Object subset) throws Exception {
244 double[] first = new double[2];
245 double[] last = new double[2];
246 int[] length = new int[2];
247
248 HashMap<String, double[]> domainSubset = new HashMap<String, double[]>();
249 domainSubset.put(trackDim_name, (double[]) ((HashMap)subset).get(trackDim_name));
250 domainSubset.put(vertDim_name, (double[]) ((HashMap)subset).get(vertDim_name));
251
252 for (int kk=0; kk<2; kk++) {
253 RealType rtype = domainRealTypes[kk];
254 double[] coords = (double[]) ((HashMap)subset).get(rtype.getName());
255 first[kk] = coords[0];
256 last[kk] = coords[1];
257 length[kk] = (int) ((last[kk] - first[kk])/coords[2] + 1);
258 last[kk] = first[kk]+coords[2]*(length[kk]-1);
259 }
260 Linear2DSet domainSet = new Linear2DSet(first[0], last[0], length[0], first[1], last[1], length[1]);
261 final Linear1DSet[] lin1DSet_s = new Linear1DSet[] {domainSet.getLinear1DComponent(0),
262 domainSet.getLinear1DComponent(1)};
263
264 float[] new_altitudes = new float[length[vert_tup_idx]];
265 float[] new_times = new float[length[track_tup_idx]];
266
267 int track_start = (int) first[track_tup_idx];
268 int vert_start = (int) first[vert_tup_idx];
269 int vert_skip = (int) ((double[]) ((HashMap)subset).get(vertDim_name))[2];
270 int track_skip = (int) ((double[]) ((HashMap)subset).get(trackDim_name))[2];
271 for (int k=0; k<new_altitudes.length; k++) {
272 new_altitudes[k] = vertLocs[vert_start+(k*vert_skip)];
273 }
274
275 for (int k=0; k<new_times.length; k++) {
276 new_times[k] = trackTimes[track_start+(k*track_skip)];
277 }
278
279 final Gridded1DSet alt_set = new Gridded1DSet(vertLocType, new float[][] {new_altitudes}, new_altitudes.length);
280 final Gridded1DSet time_set = new Gridded1DSet(trackTimeType, new float[][] {new_times}, new_times.length);
281 final float vert_offset = (float) first[vert_tup_idx];
282 final float track_offset = (float) first[track_tup_idx];
283
284 RealTupleType reference = new RealTupleType(vertLocType, trackTimeType);
285
286 CoordinateSystem cs = null;
287
288 try {
289 cs = new CoordinateSystem(reference, null) {
290 public float[][] toReference(float[][] vals) throws VisADException {
291 int[] indexes = lin1DSet_s[0].valueToIndex(new float[][] {vals[0]});
292 for (int k=0; k<vals[0].length;k++) {
293 //-indexes[k] = (int) (vals[vert_tup_idx][k] - vert_offset); ?
294 }
295 float[][] alts = alt_set.indexToValue(indexes);
296
297 indexes = lin1DSet_s[1].valueToIndex(new float[][] {vals[1]});
298 for (int k=0; k<vals[0].length;k++) {
299 //-indexes[k] = (int) (vals[track_tup_idx][k] - track_offset); ?
300 }
301 float[][] times = time_set.indexToValue(indexes);
302
303 return new float[][] {alts[0], times[0]};
304 }
305 public float[][] fromReference(float[][] vals) throws VisADException {
306 int[] indexes = alt_set.valueToIndex(new float[][] {vals[vert_tup_idx]});
307 float[][] vert_coords = lin1DSet_s[vert_tup_idx].indexToValue(indexes);
308 indexes = time_set.valueToIndex(new float[][] {vals[track_tup_idx]});
309 float[][] track_coords = lin1DSet_s[track_tup_idx].indexToValue(indexes);
310 return new float[][] {vert_coords[0], track_coords[0]};
311 }
312 public double[][] toReference(double[][] vals) throws VisADException {
313 return Set.floatToDouble(toReference(Set.doubleToFloat(vals)));
314 }
315 public double[][] fromReference(double[][] vals) throws VisADException {
316 return Set.floatToDouble(fromReference(Set.doubleToFloat(vals)));
317 }
318 public boolean equals(Object obj) {
319 return true;
320 }
321 };
322 }
323 catch (Exception e) {
324 }
325
326 RealTupleType domainTupType = new RealTupleType(domainRealTypes[0], domainRealTypes[1], cs, null);
327 domainSet = new Linear2DSet(domainTupType, first[0], last[0], length[0], first[1], last[1], length[1]);
328
329 return domainSet;
330 }
331
332 public FunctionType getMathType() {
333 return null;
334 }
335
336 public RealType[] getDomainRealTypes() {
337 return domainRealTypes;
338 }
339
340 public HashMap getDefaultSubset() {
341 HashMap subset = ProfileAlongTrack.getEmptySubset();
342
343 double[] coords = (double[])subset.get("TrackDim");
344 coords[0] = 20000.0;
345 coords[1] = (TrackLen - 15000.0) - 1;
346 //-coords[2] = 30.0;
347 coords[2] = 5.0;
348 subset.put("TrackDim", coords);
349
350 coords = (double[])subset.get("VertDim");
351 coords[0] = 98.0;
352 coords[1] = (VertLen) - 1;
353 coords[2] = 2.0;
354 subset.put("VertDim", coords);
355 return subset;
356 }
357
358 public int[] getTrackRangeInsideLonLatRect(double minLat, double maxLat, double minLon, double maxLon) {
359 int nn = 100;
360 int skip = TrackLen/nn;
361 double lon;
362 double lat;
363
364 int idx = 0;
365 while (idx < TrackLen) {
366 lon = (double)trackLongitude[idx];
367 lat = (double)trackLatitude[idx];
368 if (((lon > minLon) && (lon < maxLon)) && ((lat > minLat)&&(lat < maxLat))) break;
369 idx += skip;
370 }
371 if (idx > TrackLen-1) idx = TrackLen-1;
372 if (idx == TrackLen-1) return new int[] {-1,-1};
373
374 int low_idx = idx;
375 while (low_idx > 0) {
376 lon = (double)trackLongitude[low_idx];
377 lat = (double)trackLatitude[low_idx];
378 if (((lon > minLon) && (lon < maxLon)) && ((lat > minLat)&&(lat < maxLat))) {
379 low_idx -= 1;
380 continue;
381 }
382 else {
383 break;
384 }
385 }
386
387 int hi_idx = idx;
388 while (hi_idx < TrackLen-1) {
389 lon = (double)trackLongitude[hi_idx];
390 lat = (double)trackLatitude[hi_idx];
391 if (((lon > minLon) && (lon < maxLon)) && ((lat > minLat)&&(lat < maxLat))) {
392 hi_idx += 1;
393 continue;
394 }
395 else {
396 break;
397 }
398 }
399 return new int[] {low_idx, hi_idx};
400 }
401
402 public HashMap getSubsetFromLonLatRect(HashMap subset, double minLat, double maxLat, double minLon, double maxLon) {
403 double[] coords = (double[])subset.get("TrackDim");
404 int[] idxs = getTrackRangeInsideLonLatRect(minLat, maxLat, minLon, maxLon);
405 coords[0] = (double) idxs[0];
406 coords[1] = (double) idxs[1];
407 if ((coords[0] == -1) || (coords[1] == -1)) return null;
408 return subset;
409 }
410
411 public HashMap getSubsetFromLonLatRect(HashMap subset, double minLat, double maxLat, double minLon, double maxLon,
412 int xStride, int yStride, int zStride) {
413
414 double[] coords = (double[])subset.get("TrackDim");
415 int[] idxs = getTrackRangeInsideLonLatRect(minLat, maxLat, minLon, maxLon);
416 coords[0] = (double) idxs[0];
417 coords[1] = (double) idxs[1];
418 if ((coords[0] == -1) || (coords[1] == -1)) return null;
419 if (xStride > 0) {
420 coords[2] = xStride;
421 }
422
423 coords = (double[])subset.get("VertDim");
424 if (yStride > 0) {
425 coords[2] = yStride;
426 }
427 return subset;
428 }
429
430 public HashMap getSubsetFromLonLatRect(double minLat, double maxLat, double minLon, double maxLon) {
431 return getSubsetFromLonLatRect(getDefaultSubset(), minLat, maxLat, minLon, maxLon);
432 }
433
434 public HashMap getSubsetFromLonLatRect(double minLat, double maxLat, double minLon, double maxLon,
435 int xStride, int yStride, int zStride) {
436
437 return getSubsetFromLonLatRect(getDefaultSubset(), minLat, maxLat, minLon, maxLon, xStride, yStride, zStride);
438 }
439
440 public abstract float[] getVertBinAltitude() throws Exception;
441 public abstract float[] getTrackTimes() throws Exception;
442 public abstract RealType makeVertLocType() throws Exception;
443 public abstract RealType makeTrackTimeType() throws Exception;
444 public abstract float[] getTrackLongitude() throws Exception;
445 public abstract float[] getTrackLatitude() throws Exception;
446
447 public static float[] medianFilter(float[] A, int lenx, int leny, int window_lenx, int window_leny)
448 throws VisADException {
449 float[] result = new float[A.length];
450 float[] window = new float[window_lenx*window_leny];
451 float[] new_window = new float[window_lenx*window_leny];
452 int[] sort_indexes = new int[window_lenx*window_leny];
453
454 int a_idx;
455 int w_idx;
456
457 int w_lenx = window_lenx/2;
458 int w_leny = window_leny/2;
459
460 int lo;
461 int hi;
462 int ww_jj;
463 int ww_ii;
464 int cnt;
465
466 for (int j=0; j<leny; j++) {
467 for (int i=0; i<lenx; i++) {
468 a_idx = j*lenx + i;
469 cnt = 0;
470 for (int w_j=-w_leny; w_j<w_leny; w_j++) {
471 for (int w_i=-w_lenx; w_i<w_lenx; w_i++) {
472 ww_jj = w_j + j;
473 ww_ii = w_i + i;
474 w_idx = (w_j+w_leny)*window_lenx + (w_i+w_lenx);
475 if ((ww_jj >= 0) && (ww_ii >=0) && (ww_jj < leny) && (ww_ii < lenx)) {
476 window[cnt] = A[ww_jj*lenx+ww_ii];
477 cnt++;
478 }
479 }
480 }
481 System.arraycopy(window, 0, new_window, 0, cnt);
482 //-sort_indexes = QuickSort.sort(new_window, sort_indexes);
483 sort_indexes = QuickSort.sort(new_window);
484 result[a_idx] = new_window[cnt/2];
485 }
486 }
487 return result;
488 }
489
490 }