1

I am trying to collocate 2 datasets of different resolution using haversine's formula in python.

This loop taking too much time to compute the distance using the haversine formula. In order to reduce the unwanted collocating 2 very distant data points, I gave a condition to skip to next grid point

for i in range(len(lati1)):
for j in range (len(long1)):
    for x in range(len(lati2)):
        for y in range(len(long2)):

            lat1 = radians(lati1[i])
            lon1 = radians(long1[j])
            lat2 = radians(lati2[x])
            lon2 = radians(long2[y])



       temp[0][i][j]
            dlat = lat2 - lat1
            dlon = lon2 - lon1
            a = sin(dlat/2)**2 + cos(lat1) * cos(lat2) * sin(dlon / 2)**2
            c = 2 * atan2(sqrt(a), sqrt(1 - a))
            distance = R * c 
            print distance
            if ( distance <= 3.5):
                print distance, lati1[i], long1[j],lati2[x], long2[y], temp[i][j], preci[x][y]
                writer.writelines("\r")
                writer.writelines( (' ').join([distance, lati1[i], long1[j], temp[i][j], preci[x][y]]))
            elif (distance < 5.0):
                break
    break
j=j+1

after running this program for a long time, it is showing like this

572
572
3.2412380942 0.13999939 60.100006 0.125 60.125 287.77 0.0
Traceback (most recent call last):

  File "<ipython-input-36-ac8328bf9abd>", line 40, in <module>
    writer.writelines( (' ').join([distance, btlat[i], btlon[j], btlat[0,i,j], predata[x,y]]))

  File "/home/krishna/.local/lib/python2.7/site-packages/numpy/ma/core.py", line 3174, in __getitem__
    dout = self.data[indx]

IndexError: too many indices for array

~ I can change the array indices but is there any better solutions? to collocate?

Some could you help me to resolve this issue? Thanks

Krishnaap
  • 297
  • 3
  • 18
  • 1
    Within your `for y in range(len(prelon)):` you have a `y = y + 1` (and possibly something similar for `j`, but that is difficult to say since the indentation seems to be incorrect), what is that supposed to do? – Bart Jul 31 '19 at 13:53
  • And the traceback has both `btlat[i]` and `btlat[0,i,j]`; how many dimensions does `btlat` have? – Bart Jul 31 '19 at 13:57
  • Yes, the y=y+1, it is supposed to increment in the index. So is it possible to do that without writing it? You mean btdata? it has 3 dimensions(1, 501, 572), even though time dimension is 1 only. – Krishnaap Jul 31 '19 at 14:04
  • I mean after collocating one data point in first dataset (temp) with preci, the program will check collocation with all the data points of preci. So to skip this unwanted calculation I want to skip to next temp grid point in the next longitude so how to jump to j=j+1 ? – Krishnaap Jul 31 '19 at 14:11
  • Can you check your code/traceback to make sure it is correct (what is that `btdata[0][i][j]` doing there in the middle of the `for y in range()` loop?) and correctly indented (the `else` in the last part of the script...)? – Bart Jul 31 '19 at 14:20
  • Basically, I am collocating every data points of temp data with precipitation data. That's why I used the for loops like that – Krishnaap Jul 31 '19 at 14:27
  • Is there any way to go back to specific for loop? or skip to a specific for loop? if this line agrees elif (distance < 5.0): y =y + 1 go to the for loop for j in range (len(btlon)): – Krishnaap Jul 31 '19 at 14:30
  • 1
    Breaking out of multple loops is usually not the pretiest, but there are some solutions e.g. here: https://stackoverflow.com/questions/189645/how-to-break-out-of-multiple-loops-in-python – Bart Jul 31 '19 at 14:42
  • Thanks bart, let me check that option. – Krishnaap Jul 31 '19 at 14:48
  • I skipped to the second for loop but did not get any outputs elif (distance < 5.0): break break j=j+1 – Krishnaap Jul 31 '19 at 15:11
  • skipping from 2 for loops and giving increment j=j+1 is not working, since for y in range taking from consecutive blon values – Krishnaap Jul 31 '19 at 15:20
  • I m not getting any outputs, though the program runs for 20-40 min. – Krishnaap Jul 31 '19 at 16:05
  • It would be helpful if you include the updated code in you question. Or even more helpful if you create a [minimal working example](https://stackoverflow.com/help/minimal-reproducible-example): ditch the NetCDF files, replace the lat/lons from both files with some dummy lat/lon arrays, that way everyone can copy/paste/run/debug your code and suggest improvements. – Bart Jul 31 '19 at 16:52

0 Answers0