Source code

001/*
002 * This file is part of McIDAS-V
003 *
004 * Copyright 2007-2024
005 * Space Science and Engineering Center (SSEC)
006 * University of Wisconsin - Madison
007 * 1225 W. Dayton Street, Madison, WI 53706, USA
008 * https://www.ssec.wisc.edu/mcidas/
009 * 
010 * All Rights Reserved
011 * 
012 * McIDAS-V is built on Unidata's IDV and SSEC's VisAD libraries, and
013 * some McIDAS-V source code is based on IDV and VisAD source code.  
014 * 
015 * McIDAS-V is free software; you can redistribute it and/or modify
016 * it under the terms of the GNU Lesser Public License as published by
017 * the Free Software Foundation; either version 3 of the License, or
018 * (at your option) any later version.
019 * 
020 * McIDAS-V is distributed in the hope that it will be useful,
021 * but WITHOUT ANY WARRANTY; without even the implied warranty of
022 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
023 * GNU Lesser Public License for more details.
024 * 
025 * You should have received a copy of the GNU Lesser Public License
026 * along with this program.  If not, see https://www.gnu.org/licenses/.
027 */
028
029package edu.wisc.ssec.mcidasv.adt;
030
031import java.io.IOException;
032import java.lang.Math;
033import java.lang.String;
034import java.rmi.RemoteException;
035import java.util.List;
036
037import org.slf4j.Logger;
038import org.slf4j.LoggerFactory;
039
040import ucar.unidata.data.grid.GridUtil;
041import ucar.unidata.util.IOUtil;
042import ucar.unidata.util.StringUtil;
043import visad.FlatField;
044import visad.VisADException;
045
046public class ReadIRImage {
047
048    private static final Logger logger = LoggerFactory.getLogger(ReadIRImage.class);
049
050    public ReadIRImage() {
051    }
052
053    public static class ImgCoeffs {
054        String sat_id;
055        int sat_num;
056        int chan;
057        int det;
058        float scal_m;
059        float scal_b;
060        float side;
061        float conv_n;
062        float conv_a;
063        float conv_b;
064        float conv_g;
065
066        public ImgCoeffs(List<String> toks) {
067            sat_id = toks.get(0);
068            sat_num = Integer.parseInt(toks.get(1));
069            chan = Integer.parseInt(toks.get(2));
070            det = Integer.parseInt(toks.get(3));
071            scal_m = Float.parseFloat(toks.get(4));
072            scal_b = Float.parseFloat(toks.get(5));
073            side = Float.parseFloat(toks.get(6));
074            conv_n = Float.parseFloat(toks.get(7));
075            conv_a = Float.parseFloat(toks.get(8));
076            conv_b = Float.parseFloat(toks.get(9));
077            conv_g = Float.parseFloat(toks.get(10));
078        }
079
080        public String getSat_id() {
081            return sat_id;
082        }
083    }
084
085    public static void ReadIRDataFile(FlatField satgrid, float cenlat, float cenlon,
086            int SatelliteID, int satChannel, boolean isTemperature) throws IOException,
087            VisADException, RemoteException {
088        // Retrieve temperatures from image. This to be done in IDV
089        GridUtil.Grid2D g2d = GridUtil.makeGrid2D(satgrid);
090        float[][] temps = null;
091        float[][][] satimage = null;
092        int numx = 200;
093        int numy = 200;
094        float[][] LocalLatitude = new float[200][200];
095        float[][] LocalLongitude = new float[200][200];
096        float[][] LocalTemperature = new float[200][200];
097
098        // now spatial subset numx by numy
099        GridUtil.Grid2D g2d1 = spatialSubset(g2d, cenlat, cenlon, numx, numy);
100
101        satimage = g2d1.getvalues();
102        float[][] temp0 = satimage[0];
103
104        if (isTemperature) {
105            temps = temp0;
106        } else {
107            temps = im_gvtota(numx, numy, temp0, SatelliteID, satChannel);
108        }
109
110        Data.IRData_NumberRows = numy;
111        Data.IRData_NumberColumns = numx;
112        Data.IRData_CenterLatitude = cenlat;
113        Data.IRData_CenterLongitude = cenlon;
114
115        LocalTemperature = temps;
116        LocalLatitude = g2d1.getlats();
117        LocalLongitude = g2d1.getlons();
118
119        for (int XInc = 0; XInc < numx; XInc++) {
120            for (int YInc = 0; YInc < numy; YInc++) {
121                // must flip x/y to y/x for ADT automated routines
122                Data.IRData_Latitude[YInc][XInc] = LocalLatitude[XInc][YInc];
123                Data.IRData_Longitude[YInc][XInc] = LocalLongitude[XInc][YInc];
124                Data.IRData_Temperature[YInc][XInc] = LocalTemperature[XInc][YInc];
125            }
126        }
127        int CenterXPos = Data.IRData_NumberColumns / 2;
128        int CenterYPos = Data.IRData_NumberRows / 2;
129
130        double LocalValue[] = Functions.distance_angle(
131                Data.IRData_Latitude[CenterYPos][CenterXPos],
132                Data.IRData_Longitude[CenterYPos][CenterXPos],
133                Data.IRData_Latitude[CenterYPos + 1][CenterXPos],
134                Data.IRData_Longitude[CenterYPos][CenterXPos], 1);
135
136        Data.IRData_ImageResolution = LocalValue[0];
137
138        History.IRCurrentRecord.date = Data.IRData_JulianDate;
139        History.IRCurrentRecord.time = Data.IRData_HHMMSSTime;
140        History.IRCurrentRecord.latitude = Data.IRData_CenterLatitude;
141        History.IRCurrentRecord.longitude = Data.IRData_CenterLongitude;
142        History.IRCurrentRecord.sattype = SatelliteID;
143
144        int RetVal[] = Functions.adt_oceanbasin(Data.IRData_CenterLatitude,
145                Data.IRData_CenterLongitude);
146        Env.DomainID = RetVal[1];
147        // System.out.printf("lat=%f lon=%f domainID=%d\n",Data.IRData_CenterLatitude,Data.IRData_CenterLongitude,Env.DomainID);
148    }
149
150    private static GridUtil.Grid2D spatialSubset(GridUtil.Grid2D g2d, float cenlat, float cenlon,
151            int numx, int numy) {
152        float[][] lats = g2d.getlats();
153        float[][] lons = g2d.getlons();
154        float[][][] values = g2d.getvalues();
155        float[][] slats = new float[numx][numy];
156        float[][] slons = new float[numx][numy];
157        float[][][] svalues = new float[1][numx][numy];
158
159        int ly = lats[0].length;
160        int ly0 = ly / 2;
161        int lx = lats.length;
162        logger.debug("lenx: {}, leny: {}", lx, ly);
163        int lx0 = lx / 2;
164        int ii = numx / 2, jj = numy / 2;
165
166        for (int j = 0; j < ly - 1; j++) {
167            if (Float.isNaN(lats[lx0][j])) {
168                continue;
169            }
170            if ((lats[lx0][j] > cenlat) && (lats[lx0][j + 1] < cenlat)) {
171                jj = j;
172            }
173        }
174        for (int i = 0; i < lx - 1; i++) {
175            if (Float.isNaN(lons[i][ly0])) {
176                continue;
177            }
178            if ((lons[i][ly0] < cenlon) && (lons[i + 1][ly0] > cenlon)) {
179                ii = i;
180            }
181        }
182        int startx = ii - (numx / 2 - 1);
183        int starty = jj - (numy / 2 - 1);
184        logger.debug("startx: {}, starty: {}", startx, starty);
185        logger.debug("numx: {}, numy: {}", numx, numy);
186
187        if (startx < 0) {
188            startx = 0;
189        }
190        if (starty < 0) {
191            starty = 0;
192        }
193        for (int i = 0; i < numx; i++) {
194            for (int j = 0; j < numy; j++) {
195                try {
196                    slats[i][j] = lats[i + startx][j + starty];
197                } catch (ArrayIndexOutOfBoundsException aioobe) {
198                    slats[i][j] = Float.NaN;
199                }
200                try {
201                    slons[i][j] = lons[i + startx][j + starty];
202                } catch (ArrayIndexOutOfBoundsException aioobe) {
203                    slats[i][j] = Float.NaN;
204                }
205                try {
206                    svalues[0][i][j] = values[0][i + startx][j + starty];
207                } catch (ArrayIndexOutOfBoundsException aioobe) {
208                    slats[i][j] = Float.NaN;
209                }
210            }
211        }
212        return new GridUtil.Grid2D(slats, slons, svalues);
213    }
214
215    /**
216     * im_gvtota
217     *
218     * This subroutine converts GVAR counts to actual temperatures based on the
219     * current image set in IM_SIMG.
220     *
221     * im_gvtota ( int *nvals, unsigned int *gv, float *ta, int *iret )
222     *
223     * Input parameters: *nvals int Number of values to convert *gv int Array of
224     * GVAR count values
225     *
226     * Output parameters: *ta float Array of actual temperatures *iret int
227     * Return value = -1 - could not open table = -2 - could not find match
228     *
229     *
230     * Log: D.W.Plummer/NCEP 02/03 D.W.Plummer/NCEP 06/03 Add coeff G for 2nd
231     * order poly conv T. Piper/SAIC 07/06 Added tmpdbl to eliminate warning
232     */
233    private static float[][] im_gvtota(int nx, int ny, float[][] gv, int imsorc, int imtype) {
234        double c1 = 1.191066E-5;
235        double c2 = 1.438833;
236
237        int ii;
238        int ip;
239        int chan;
240        int found;
241        double Rad;
242        double Teff;
243        double tmpdbl;
244        float[][] ta = new float[nx][ny];
245        String fp = "/ucar/unidata/data/storm/ImgCoeffs.tbl";
246
247        for (ii = 0; ii < nx; ii++) {
248            for (int jj = 0; jj < ny; jj++) {
249                ta[ii][jj] = Float.NaN;
250            }
251        }
252
253        // Read in coefficient table if necessary.
254        String s = null;
255        try {
256            s = IOUtil.readContents(fp);
257        } catch (Exception re) {
258            logger.warn("Failed to read coefficient table", re);
259        }
260
261        int i = 0;
262        ImgCoeffs[] ImageConvInfo = new ImgCoeffs[50];
263        for (String line : StringUtil.split(s, "\n", true, true)) {
264            if (line.startsWith("!")) {
265                continue;
266            }
267            List<String> stoks = StringUtil.split(line, " ", true, true);
268            ImageConvInfo[i] = new ImgCoeffs(stoks);
269            i++;
270        }
271        int nImgRecs = i;
272        found = 0;
273        ii = 0;
274        while ((ii < nImgRecs) && (found == 0)) {
275            int tmp = ImageConvInfo[ii].chan - 1;
276            tmpdbl = (double) (tmp * tmp);
277            chan = G_NINT(tmpdbl);
278            if ((imsorc == ImageConvInfo[ii].sat_num) && (imtype == chan)) {
279                found = 1;
280            } else {
281                ii++;
282            }
283        }
284
285        if (found == 0) {
286            return null;
287        } else {
288            ip = ii;
289            for (ii = 0; ii < nx; ii++) {
290                for (int jj = 0; jj < ny; jj++) {
291                    // Convert GVAR count (gv) to Scene Radiance
292                    Rad = ((double) gv[ii][jj] - ImageConvInfo[ip].scal_b) /
293                    // -------------------------------------
294                            ImageConvInfo[ip].scal_m;
295
296                    Rad = Math.max(Rad, 0.0);
297
298                    // Convert Scene Radiance to Effective Temperature
299                    Teff = (c2 * ImageConvInfo[ip].conv_n) /
300                    /*
301                     * -------------------------------------------------- -----
302                     */
303                    (Math.log(1.0 + (c1 * Math.pow(ImageConvInfo[ip].conv_n, 3.0)) / Rad));
304                    // Convert Effective Temperature to Temperature
305                    ta[ii][jj] = (float) (ImageConvInfo[ip].conv_a + ImageConvInfo[ip].conv_b
306                            * Teff + ImageConvInfo[ip].conv_g * Teff * Teff);
307                }
308            }
309        }
310        return ta;
311    }
312
313    public static int G_NINT(double x) {
314        return (((x) < 0.0F) ? ((((x) - (float) ((int) (x))) <= -.5f) ? (int) ((x) - .5f)
315                : (int) (x)) : ((((x) - (float) ((int) (x))) >= .5f) ? (int) ((x) + .5f)
316                : (int) (x)));
317    }
318}